text:
@muladd function perform_step!(integrator,cache::DRI1ConstantCache,f=integrator.f)
  @unpack a021,a031,a032,a121,a131,b021,b031,b121,b131,b221,b222,b223,b231,b232,b233,α1,α2,α3,c02,c03,c12,c13,beta11,beta12,beta13,beta22,beta23,beta31,beta32,beta33,beta42,beta43,NORMAL_ONESIX_QUANTILE = cache
  @unpack t,dt,uprev,u,W,p = integrator

  # define three-point distributed random variables
  dW_scaled = W.dW / sqrt(dt)
  _dW = map(x -> calc_threepoint_random(integrator, NORMAL_ONESIX_QUANTILE, x), dW_scaled)
  chi1 = map(x -> (x^2-dt)/2, _dW) # diagonal of Ihat2
  if !(typeof(W.dW) <: Number)
    m = length(W.dW)
    # define two-point distributed random variables
    _dZ = map(x -> calc_twopoint_random(integrator, x),  W.dZ)
    Ihat2 = zeros(eltype(W.dZ), m, m) # I^_(k,l)
    for k = 1:m
      for l = 1:m
        if k<l
          Ihat2[k, l] = (_dW[k]*_dW[l]-integrator.sqdt*_dZ[k])/2
        elseif l<k
          Ihat2[k, l] = (_dW[k]*_dW[l]+integrator.sqdt*_dZ[l])/2
        else k==l
          Ihat2[k, k] = chi1[k]
        end
      end
    end
  end

  # compute stage values
  k1 = integrator.f(uprev,p,t)
  g1 = integrator.g(uprev,p,t)

  # H_i^(0), stage 1
  # H01 = uprev
  # H_i^(0), stage 2
  if !is_diagonal_noise(integrator.sol.prob) || typeof(W.dW) <: Number
    H02 = uprev + a021*k1*dt + b021*g1*_dW
  else
    H02 = uprev + a021*k1*dt + b021*g1.*_dW
  end

  # H_i^(0), stage 3 (requires H_i^(k) stage 2 in general)
  k2 = integrator.f(H02,p,t+c02*dt)
  H03 = uprev + a032*k2*dt + a031*k1*dt
  if !is_diagonal_noise(integrator.sol.prob) || typeof(W.dW) <: Number
    H03 += b031*g1*_dW
  else
    H03 += b031*g1.*_dW
  end

  k3 = integrator.f(H03,p,t+c03*dt)

  # H_i^(k), stage 1
  # H11 = uprev
  # H_i^(k), stage 2
  if typeof(W.dW) <: Number
    H12 = uprev + a121*k1*dt + b121*g1*integrator.sqdt
  elseif is_diagonal_noise(integrator.sol.prob)
    H12 = Vector{typeof(uprev)}[uprev .+ a121*k1*dt .+ b121*g1[k]*integrator.sqdt for k=1:m]
  else
    H12 = Vector{typeof(uprev)}[uprev .+ a121*k1*dt .+ b121*g1[:,k]*integrator.sqdt for k=1:m]
  end

  # # H_i^(k), stage 3
  if typeof(W.dW) <: Number
    H13 = uprev + a131*k1*dt + b131*g1*integrator.sqdt
  elseif is_diagonal_noise(integrator.sol.prob)
    H13 = Vector{typeof(uprev)}[uprev .+ a131*k1*dt .+ b131*g1[k]*integrator.sqdt for k=1:m]
  else
    H13 = Vector{typeof(uprev)}[uprev .+ a131*k1*dt .+ b131*g1[:,k]*integrator.sqdt for k=1:m]
  end

  # H^_i^(k), stage 1
  # H21 = uprev
  # H^_i^(k), stage 2

  if typeof(W.dW) <: Number
    g2 = integrator.g(H12,p,t+c12*dt)
    g3 = integrator.g(H13,p,t+c13*dt)
    # for m=1:  H22 = uprev
  else
    g2 = [integrator.g(H12[k],p,t+c12*dt) for k=1:m]
    g3 = [integrator.g(H13[k],p,t+c13*dt) for k=1:m]
    H22 = [uprev for k=1:m]
    # add inbounds for speed if working properly
    for k=1:m
      for l=1:m
        if k == l
          continue
        end
        if is_diagonal_noise(integrator.sol.prob)
          @.. H22[k] += (b221*g1[l]+b222*g2[l][l]+b223*g3[l][l])*Ihat2[k,l]/integrator.sqdt
        else
          H22[k] += (b221*g1[:,l]+b222*g2[l][:,l]+b223*g3[l][:,l])*Ihat2[k,l]/integrator.sqdt
        end
      end
    end
  end

  # H^_i^(k), stage 3
  # for m=1: H23 = uprev
  if !(typeof(W.dW) <: Number)
    H23 = [uprev for k=1:m]
    # add inbounds for speed if working properly
    for k=1:m
      for l=1:m
        if k == l
          continue
        end
        if is_diagonal_noise(integrator.sol.prob)
          @.. H23[k] += (b231*g1[l]+b232*g2[l][l]+b233*g3[l][l])*Ihat2[k,l]/integrator.sqdt
        else
          H23[k] += (b231*g1[:,l]+b232*g2[l][:,l]+b233*g3[l][:,l])*Ihat2[k,l]/integrator.sqdt
        end
      end
    end
  end

  # add stages together Eq. (3)
  u = uprev + α1*k1*dt + α2*k2*dt + α3*k3*dt

  # add noise
  if typeof(W.dW) <: Number
    # lines 2 and 3
    u += g1*(_dW*beta11) + g2*(_dW*beta12+chi1*beta22/integrator.sqdt) + g3*(_dW*beta13+chi1*beta23/integrator.sqdt)
    # lines 4 and 5 are zero by construction
    # u += g1*(_dW*(beta31+beta32+beta33)+chi1*integrator.sqdt*(beta42+beta43))
  else
      if is_diagonal_noise(integrator.sol.prob)
        u += g1.*_dW*beta11
        for k=1:m
          u[k] += g2[k][k]*_dW[k]*beta12+g3[k][k]*_dW[k]*beta13
          u[k] += g2[k][k]*chi1[k]*beta22/integrator.sqdt + g3[k][k]*chi1[k]*beta23/integrator.sqdt
          tmpg = integrator.g(H22[k],p,t)
          u[k] = u[k] + g1[k]*_dW[k]*beta31 + tmpg[k]*(_dW[k]*beta32 + integrator.sqdt*beta42)
          tmpg = integrator.g(H23[k],p,t)
          u[k] = u[k] + tmpg[k]*(_dW[k]*beta33 + integrator.sqdt*beta43)
        end
      else
        # non-diag noise
        for k=1:m
          g1k = @view g1[:,k]
          g2k = @view g2[k][:,k]
          g3k = @view g3[k][:,k]
          @.. u = u + g1k*_dW[k]*(beta11+beta31)
          @.. u = u + g2k*_dW[k]*beta12 + g3k*_dW[k]*beta13
          @.. u = u + g2k*chi1[k]*beta22/integrator.sqdt + g3k*chi1[k]*beta23/integrator.sqdt
          tmpg = integrator.g(H22[k],p,t)
          @.. u = u + tmpgk*_dW[k]*beta32 + tmpgk*integrator.sqdt*beta42
          tmpg = integrator.g(H23[k],p,t)
          @.. u = u + tmpgk*_dW[k]*beta33 + tmpgk*integrator.sqdt*beta43
        end
      end
  end

  # Test EM scheme:
  # integrator.u = uprev + 1*k1*dt + 1*g1*_dW

  integrator.u = u

end

@muladd function perform_step!(integrator,cache::DRI1Cache,f=integrator.f)
  @unpack t,dt,uprev,u,W,p = integrator
  @unpack _dW,_dZ,chi1,Ihat2,tab,g1,g2,g3,k1,k2,k3,H02,H03,H12,H13,H22,H23,tmp1,tmpg = cache
  @unpack a021,a031,a032,a121,a131,b021,b031,b121,b131,b221,b222,b223,b231,b232,b233,α1,α2,α3,c02,c03,c12,c13,beta11,beta12,beta13,beta22,beta23,beta31,beta32,beta33,beta42,beta43,NORMAL_ONESIX_QUANTILE = cache.tab

  m = length(W.dW)

  if typeof(W.dW) <: Union{SArray,Number}
    # tbd
    _dW = map(x -> calc_threepoint_random(integrator, NORMAL_ONESIX_QUANTILE, x), W.dW / sqrt(dt))
  else
    # define three-point distributed random variables
    @.. chi1 = W.dW / sqrt(dt)
    calc_threepoint_random!(_dW, integrator, NORMAL_ONESIX_QUANTILE, chi1)
    map!(x -> (x^2-dt)/2, chi1, _dW)

    # define two-point distributed random variables
    calc_twopoint_random!(_dZ,integrator,W.dZ)

    for k = 1:m
      for l = 1:m
        if k<l
          Ihat2[k, l] = (_dW[k]*_dW[l]-integrator.sqdt*_dZ[k])/2
        elseif l<k
          Ihat2[k, l] = (_dW[k]*_dW[l]+integrator.sqdt*_dZ[l])/2
        else k==l
          Ihat2[k, k] = chi1[k]
        end
      end
    end
  end

  # compute stage values
  integrator.f(k1,uprev,p,t)
  integrator.g(g1,uprev,p,t)

  # H_i^(0), stage 1
  # H01 = uprev
  # H_i^(0), stage 2
  if is_diagonal_noise(integrator.sol.prob)
    @.. H02 = uprev + a021*k1*dt + b021*g1*_dW
  else
    mul!(tmp1,g1,_dW)
    @.. H02 = uprev + dt*a021*k1 + b021*tmp1
  end

  integrator.f(k2,H02,p,t+c02*dt)
  if is_diagonal_noise(integrator.sol.prob)
    @.. H03 = uprev + a032*k2*dt + a031*k1*dt + b031*g1*_dW
  else
    @.. H03 = uprev + a032*k2*dt + a031*k1*dt + b031*tmp1
  end

  integrator.f(k3,H03,p,t+c03*dt)

  # H_i^(k), stages
  # H11 = uprev
  for k=1:m
    if is_diagonal_noise(integrator.sol.prob)
      @.. H12[k] = uprev + a121*k1*dt + b121*g1[k]*integrator.sqdt
      @.. H13[k] = uprev + a131*k1*dt + b131*g1[k]*integrator.sqdt
    else
      g1k = @view g1[:,k]
      @.. H12[k] = uprev + a121*k1*dt + b121*g1k*integrator.sqdt
      @.. H13[k] = uprev + a131*k1*dt + b131*g1k*integrator.sqdt
    end
  end
  # H^_i^(k), stages

  for k=1:m
    integrator.g(g2[k],H12[k],p,t+c12*dt)
    integrator.g(g3[k],H13[k],p,t+c13*dt)

    H22[k] = uprev
    H23[k] = uprev
  end

  for k=1:m
    for l=1:m
      if k == l
        continue
      end
      if is_diagonal_noise(integrator.sol.prob)
        @.. H22[k] = H22[k]+(b221*g1[l]+b222*g2[l][l]+b223*g3[l][l])*Ihat2[k,l]/integrator.sqdt
        @.. H23[k] = H23[k]+(b231*g1[l]+b232*g2[l][l]+b233*g3[l][l])*Ihat2[k,l]/integrator.sqdt
      else
        g1l = @view g1[:,l]
        g2l = @view g2[l][:,l]
        g3l = @view g3[l][:,l]
        @.. H22[k] = H22[k]+(b221*g1l+b222*g2l+b223*g3l)*Ihat2[k,l]/integrator.sqdt
        @.. H23[k] = H23[k]+(b231*g1l+b232*g2l+b233*g3l)*Ihat2[k,l]/integrator.sqdt
      end
    end
  end

  # add stages together Eq. (3)
  @.. u = uprev + α1*k1*dt + α2*k2*dt + α3*k3*dt

  # add noise
  if is_diagonal_noise(integrator.sol.prob)
      @.. u = u + g1*_dW*beta11
      for k=1:m
        u[k] = u[k] + g2[k][k]*_dW[k]*beta12 + g3[k][k]*_dW[k]*beta13
        u[k] = u[k]  + g2[k][k]*chi1[k]*beta22/integrator.sqdt + g3[k][k]*chi1[k]*beta23/integrator.sqdt
        if m > 1
          integrator.g(tmpg,H22[k],p,t)
          u[k] = u[k] + g1[k]*_dW[k]*beta31 + tmpg[k]*(_dW[k]*beta32 + integrator.sqdt*beta42)
          integrator.g(tmpg,H23[k],p,t)
          u[k] = u[k] + tmpg[k]*(_dW[k]*beta33 + integrator.sqdt*beta43)
        end
      end
  else
    for k=1:m
      g1k = @view g1[:,k]
      g2k = @view g2[k][:,k]
      g3k = @view g3[k][:,k]
      tmpgk = @view tmpg[:,k]
      @.. u = u + g1k*_dW[k]*(beta11+beta31)
      @.. u = u + g2k*_dW[k]*beta12 + g3k*_dW[k]*beta13
      @.. u = u + g2k*chi1[k]*beta22/integrator.sqdt + g3k*chi1[k]*beta23/integrator.sqdt
      integrator.g(tmpg,H22[k],p,t)
      @.. u = u + tmpgk*_dW[k]*beta32 + tmpgk*integrator.sqdt*beta42
      integrator.g(tmpg,H23[k],p,t)
      @.. u = u + tmpgk*_dW[k]*beta33 + tmpgk*integrator.sqdt*beta43
    end
  end

end


text:
module Walk1Ds

using POMDPs, Distributions, AdaptiveStressTesting2

export
    Walk1DParams,
    Walk1D,
    Walk1DState,
    Walk1DAction

struct Walk1DParams
    startx::Float64
    threshx::Float64
    t_max::Int
end
Walk1DParams() = Walk1DParams(1.0, 10.0, 20)

struct Walk1DState
    t::Int
    x::Float64
end
mutable struct Walk1D <: MDP{Walk1DState,Void}
    p::Walk1DParams
    distr::Distribution
    t::Int
    x::Float64
end
Walk1D(p::Walk1DParams, distr::Distribution) = Walk1D(p, distr, 0, 0.0)
Walk1D() = Walk1D(Walk1DParams(), Normal(0.0, 1.0))

POMDPs.actions(mdp::Walk1D) = [nothing]
function POMDPs.initial_state(mdp::Walk1D, rng::AbstractRNG)
    mdp.t, mdp.x = 0, mdp.p.startx
    Walk1DState(mdp.t, mdp.x)
end
function POMDPs.generate_s(mdp::Walk1D, s::Walk1DState, a::Void, rng::AbstractRNG)
    mdp.t += 1
    mdp.x += rand(rng, mdp.distr)
    Walk1DState(mdp.t, mdp.x)
end
POMDPs.isterminal(mdp::Walk1D, s::Walk1DState) = isevent(mdp, s) || mdp.t >= mdp.p.t_max
POMDPs.reward(mdp::Walk1D, s::Walk1DState, a::Void, sp::Walk1DState) = 0.0

AdaptiveStressTesting2.isevent(mdp::Walk1D, s::Walk1DState) = mdp.x > mdp.p.threshx
AdaptiveStressTesting2.state_distance(mdp::Walk1D, s1::Walk1DState, s2::Walk1DState) = abs(s1.x-s2.x)
function AdaptiveStressTesting2.miss_distance(mdp::Walk1D, s::Walk1DState)
    max(mdp.p.threshx-abs(mdp.x), 0.0)
end

Base.hash(mdp::Walk1D) = hash(mdp.t, hash(mdp.x))

function Distributions.pdf(mdp::Walk1D, s::Walk1DState, a::Void, sp::Walk1DState)
    pdf(mdp.distr, sp.x-s.x)
end

end #module

ext: jl


text:
push!(LOAD_PATH, "C:\\Users\\immanueldiener\\.julia\\dev")
using Documenter
using AutoLandmarking

makedocs(sitename="AutoLandmarking.jl")

max_stars_repo_name: imi-diener/AutoLM

max_issues_repo_name: imi-diener/AutoLM

max_forks_repo_name: imi-diener/AutoLM


text:
"""
The functions provided by ESDL are supposed to work on different types of cubes. This module defines the interface for all
Data types that
"""
module Cubes
export Axes, AbstractCubeData, getSubRange, readcubedata, AbstractCubeMem, axesCubeMem,CubeAxis, TimeAxis, TimeHAxis, QuantileAxis, VariableAxis, LonAxis, LatAxis, CountryAxis, SpatialPointAxis, caxes,
       AbstractSubCube, CubeMem, EmptyCube, YearStepRange, _read, saveCube, loadCube, RangeAxis, CategoricalAxis, axVal2Index, MSCAxis,
       getSingVal, ScaleAxis, axname, @caxis_str, rmCube, cubeproperties, findAxis, AxisDescriptor, get_descriptor, ByName, ByType, ByValue, ByFunction, getAxis,
       getOutAxis, Cube, (..), getCubeData, subsetcube, CubeMask, renameaxis!, Dataset, S3Cube, cubeinfo

"""
    AbstractCubeData{T,N}

Supertype of all cubes. `T` is the data type of the cube and `N` the number of
dimensions. Beware that an `AbstractCubeData` does not implement the `AbstractArray`
interface. However, the `ESDL` functions [`mapCube`](@ref), [`reduceCube`](@ref),
[`readcubedata`](@ref), [`plotMAP`](@ref) and [`plotXY`](@ref) will work on any subtype
of `AbstractCubeData`
"""
abstract type AbstractCubeData{T,N} end

"""
getSingVal reads a single point from the cube's data
"""
getSingVal(c::AbstractCubeData,a...)=error("getSingVal called in the wrong way with argument types $(typeof(c)), $(map(typeof,a))")

Base.eltype(::AbstractCubeData{T}) where T = T
Base.ndims(::AbstractCubeData{<:Any,N}) where N = N

"""
    readCubeData(cube::AbstractCubeData)

Given any type of `AbstractCubeData` returns a [`CubeMem`](@ref) from it.
"""
function readcubedata(x::AbstractCubeData{T,N}) where {T,N}
  s=size(x)
  aout = zeros(Union{T,Missing},s...)
  r=CartesianIndices(s)
  _read(x,aout,r)
  CubeMem(collect(CubeAxis,caxes(x)),aout,cubeproperties(x))
end

"""
This function calculates a subset of a cube's data
"""
function subsetcube end

function _read end

getsubset(x::AbstractCubeData) = x.subset === nothing ? ntuple(i->Colon(),ndims(x)) : x.subset
#"""
#Internal function to read a range from a datacube
#"""
#_read(c::AbstractCubeData,d,r::CartesianRange)=error("_read not implemented for $(typeof(c))")

"Returns the axes of a Cube"
caxes(c::AbstractCubeData)=error("Axes function not implemented for $(typeof(c))")

cubeproperties(::AbstractCubeData)=Dict{String,Any}()

"Chunks, if given"
cubechunks(c::AbstractCubeData) = (size(c,1),map(i->1,2:ndims(c))...)

"Offset of the first chunk"
chunkoffset(c::AbstractCubeData) = ntuple(i->0,ndims(c))

function iscompressed end

"Supertype of all subtypes of the original data cube"
abstract type AbstractSubCube{T,N} <: AbstractCubeData{T,N} end

"Supertype of all in-memory representations of a data cube"
abstract type AbstractCubeMem{T,N} <: AbstractCubeData{T,N} end

include("Axes.jl")
using .Axes
import .Axes: getOutAxis, getAxis
struct EmptyCube{T}<:AbstractCubeData{T,0} end
caxes(c::EmptyCube)=CubeAxis[]

"""
    CubeMem{T,N} <: AbstractCubeMem{T,N}

An in-memory data cube. It is returned by applying `mapCube` when
the output cube is small enough to fit in memory or by explicitly calling
[`readCubeData`](@ref) on any type of cube.

### Fields

* `axes` a `Vector{CubeAxis}` containing the Axes of the Cube
* `data` N-D array containing the data

"""
mutable struct CubeMem{T,N} <: AbstractCubeMem{T,N}
  axes::Vector{CubeAxis}
  data::Array{T,N}
  properties::Dict{String}
end

CubeMem(axes::Vector{CubeAxis},data) = CubeMem(axes,data,Dict{String,Any}())
Base.permutedims(c::CubeMem,p)=CubeMem(c.axes[collect(p)],permutedims(c.data,p))
caxes(c::CubeMem)=c.axes
cubeproperties(c::CubeMem)=c.properties
iscompressed(c::AbstractCubeMem)=false

Base.IndexStyle(::CubeMem)=Base.LinearFast()
function Base.setindex!(c::CubeMem,i::Integer,v)
  c.data[i] = v
end
Base.size(c::CubeMem)=size(c.data)
Base.size(c::CubeMem,i)=size(c.data,i)
Base.similar(c::CubeMem)=CubeMem(c.axes,similar(c.data))
Base.ndims(c::CubeMem{T,N}) where {T,N}=N

readcubedata(c::CubeMem)=c

function getSubRange(c::AbstractArray,i...;write::Bool=true)
  length(i)==ndims(c) || error("Wrong number of view arguments to getSubRange. Cube is: $c \n indices are $i")
  return view(c,i...)
end
getSubRange(c::Tuple{AbstractArray{T,0},AbstractArray{UInt8,0}};write::Bool=true) where {T}=c

function _subsetcube end

function subsetcube(z::CubeMem{T};kwargs...) where T
  newaxes, substuple = _subsetcube(z,collect(Any,map(Base.OneTo,size(z)));kwargs...)
  newdata = z.data[substuple...]
  !isa(newdata,AbstractArray) && (newdata = fill(newdata))
  if haskey(z.properties,"labels")
    alll = filter!(!ismissing,unique(newdata))
    z.properties["labels"] = filter(i->in(i[1],alll),z.properties["labels"])
  end
  CubeMem{T,length(newaxes)}(newaxes,newdata,z.properties)
end

"""
    gethandle(c::AbstractCubeData, [block_size])

Returns an indexable handle to the data.
"""
gethandle(c::AbstractCubeMem) = c.data
gethandle(c::CubeAxis) = collect(c.values)
gethandle(c,block_size)=gethandle(c)

import ..ESDLTools.toRange
#Generic fallback method for _read
function _read(c::CubeMem,thedata::AbstractArray,r::CartesianIndices)
  thedata .= view(c.data,r.indices...)
end

function _write(c::CubeMem,thedata::AbstractArray,r::CartesianIndices)
  cubeview = getSubRange(c.data,r.indices...)
  copyto!(cubeview,thedata)
end

#Implement getindex on AbstractCubeData objects
const IndR = Union{Integer,Colon,UnitRange}
getfirst(i::Integer,a::CubeAxis)=i
getlast(i::Integer,a::CubeAxis)=i
getfirst(i::UnitRange,a::CubeAxis)=first(i)
getlast(i::UnitRange,a::CubeAxis)=last(i)
getfirst(i::Colon,a::CubeAxis)=1
getlast(i::Colon,a::CubeAxis)=length(a)

function Base.getindex(c::AbstractCubeData,i::Integer...)
  length(i)==ndims(c) || error("You must provide $(ndims(c)) indices")
  r = CartesianIndices((first(i):first(i),Base.tail(i)...))
  aout = zeros(eltype(c),size(r))
  _read(c,aout,r)
  return aout[1]
end

function Base.getindex(c::AbstractCubeData,i::IndR...)
  length(i)==ndims(c) || error("You must provide $(ndims(c)) indices")
  ax = (caxes(c)...,)
  r = CartesianIndices(map((ii,iax)->getfirst(ii,iax):getlast(ii,iax),i,ax))
  lall = map((rr,ii)->(length(rr),!isa(ii,Integer)),r.indices,i)
  lshort = filter(ii->ii[2],collect(lall))
  newshape = map(ii->ii[1],lshort)
  aout = Array{eltype(c)}(undef,newshape...)
  if newshape == size(r)
    _read(c,aout,r)
  else
    _read(c,reshape(aout,size(r)),r)
  end
  aout
end
Base.read(d::AbstractCubeData) = getindex(d,fill(Colon(),ndims(d))...)

function formatbytes(x)
  exts=["bytes","KB","MB","GB","TB"]
  i=1
  while x>=1024
    i=i+1
    x=x/1024
  end
  return string(round(x, digits=2)," ",exts[i])
end
cubesize(c::AbstractCubeData{T}) where {T}=(sizeof(T)+1)*prod(map(length,caxes(c)))
cubesize(c::AbstractCubeData{T,0}) where {T}=sizeof(T)+1

getCubeDes(c::AbstractSubCube)="Data Cube view"
getCubeDes(::CubeAxis)="Cube axis"
getCubeDes(c::CubeMem)="In-Memory data cube"
getCubeDes(c::EmptyCube)="Empty Data Cube (placeholder)"
function Base.show(io::IO,c::AbstractCubeData)
    println(io,getCubeDes(c), " with the following dimensions")
    for a in caxes(c)
        println(io,a)
    end

    foreach(cubeproperties(c)) do p
      if p[1] in ("labels","name","units")
        println(io,p[1],": ",p[2])
      end
    end
    println(io,"Total size: ",formatbytes(cubesize(c)))
end

import ..ESDL.workdir
using NetCDF

function check_overwrite(newfolder, overwrite)
  if isdir(newfolder)
    if overwrite
      rm(newfolder, recursive=true)
    else
      error("$(newfolder) already exists, please pick another name or use `overwrite=true`")
    end
  end
end
function getsavefolder(name)
  if isempty(name)
    name = tempname()[2:end]
  end
  isabspath(name) ? name : joinpath(workdir[],name)
end

"""
    saveCube(cube,name::String)

Save a [`ZarrCube`](@ref) or [`CubeMem`](@ref) to the folder `name` in the ESDL working directory.

See also [`loadCube`](@ref)
"""
function saveCube end

# function saveCube(c::CubeMem{T},name::AbstractString; overwrite=true) where T
#   newfolder=getsavefolder(name)
#   check_overwrite(newfolder, overwrite)
#   tc=Cubes.ESDLZarr.ZArrayCube(c.axes,folder=newfolder,T=T)
#   _write(tc,c.data,CartesianIndices(c.data))
# end

import Base.Iterators: take, drop
Base.show(io::IO,a::RangeAxis)=print(io,rpad(Axes.axname(a),20," "),"Axis with ",length(a)," Elements from ",first(a.values)," to ",last(a.values))
function Base.show(io::IO,a::CategoricalAxis)
  print(io,rpad(Axes.axname(a),20," "), "Axis with ", length(a), " elements: ")
  if length(a.values)<10
    for v in a.values
      print(io,v," ")
    end
  else
    for v in take(a.values,2)
      print(io,v," ")
    end
    print(io,".. ")
    for v in drop(a.values,length(a.values)-2)
      print(io,v," ")
    end
  end
end
Base.show(io::IO,a::SpatialPointAxis)=print(io,"Spatial points axis with ",length(a.values)," points")

include("TransformedCubes.jl")
function S3Cube end
include("ZarrCubes.jl")
import .ESDLZarr: (..), Cube, getCubeData, loadCube, rmCube, CubeMask, cubeinfo, getsubinds
include("NetCDFCubes.jl")
include("Datasets.jl")
import .Datasets: Dataset, ESDLDataset
include("OBS.jl")
end


text:
using ROMEO
using Test
using NIfTI

nan_test(I1, I2) = I1[.!isnan.(I1)] ≈ I2[.!isnan.(I2)]

@testset "ROMEO.jl" begin
    include("features.jl")
    include("specialcases.jl")
    include("dsp_tests.jl")
    include("mri.jl")
    #include("timing.jl")
end


text:
recurs_sym_type(ex::Any) =
  (ex==None || typeof(ex)==Symbol || length(ex.args)==1) ? eval(ex) : Expr(ex.head, ex.args[1], recurs_sym_type(ex.args[2]))
macro c(ret_type, func, arg_types, lib)
  local _arg_types = Expr(:tuple, [recurs_sym_type(a) for a in arg_types.args]...)
  local _ret_type = recurs_sym_type(ret_type)
  local _args_in = Any[ symbol(string('a',x)) for x in 1:length(_arg_types.args) ]
  local _lib = eval(lib)
  quote
    $(esc(func))($(_args_in...)) = ccall( ($(string(func)), $(Expr(:quote, _lib)) ), $_ret_type, $_arg_types, $(_args_in...) )
  end
end

macro ctypedef(fake_t,real_t)
  real_t = recurs_sym_type(real_t)
  quote
    typealias $fake_t $real_t
  end
end

# header: /usr/local/include/cvodes/cvodes_band.h
@ctypedef CVDlsDenseJacFnB Ptr{:Void}
@ctypedef CVDlsBandJacFnB Ptr{:Void}
@c Int32 CVDlsSetDenseJacFnB (Ptr{:None},:Int32,:CVDlsDenseJacFnB) shlib
@c Int32 CVDlsSetBandJacFnB (Ptr{:None},:Int32,:CVDlsBandJacFnB) shlib
@c Int32 CVBandB (Ptr{:None},:Int32,:Clong,:Clong,:Clong) shlib

# header: /usr/local/include/cvodes/cvodes_bandpre.h
@c Int32 CVBandPrecInitB (Ptr{:None},:Int32,:Clong,:Clong,:Clong) shlib

# header: /usr/local/include/cvodes/cvodes_bbdpre.h
@ctypedef CVLocalFnB Ptr{:Void}
@ctypedef CVCommFnB Ptr{:Void}
@c Int32 CVBBDPrecInitB (Ptr{:None},:Int32,:Clong,:Clong,:Clong,:Clong,:Clong,:realtype,:CVLocalFnB,:CVCommFnB) shlib
@c Int32 CVBBDPrecReInitB (Ptr{:None},:Int32,:Clong,:Clong,:realtype) shlib

# header: /usr/local/include/cvodes/cvodes_dense.h
@c Int32 CVDenseB (Ptr{:None},:Int32,:Clong) shlib

# header: /usr/local/include/cvodes/cvodes_diag.h
@c Int32 CVDiagB (Ptr{:None},:Int32) shlib

# header: /usr/local/include/cvodes/cvodes_direct.h

# header: /usr/local/include/cvodes/cvodes.h
@ctypedef CVQuadRhsFn Ptr{:Void}
@ctypedef CVSensRhsFn Ptr{:Void}
@ctypedef CVSensRhs1Fn Ptr{:Void}
@ctypedef CVQuadSensRhsFn Ptr{:Void}
@ctypedef CVRhsFnB Ptr{:Void}
@ctypedef CVRhsFnBS Ptr{:Void}
@ctypedef CVQuadRhsFnB Ptr{:Void}
@ctypedef CVQuadRhsFnBS Ptr{:Void}
@c Int32 CVodeQuadInit (Ptr{:None},:CVQuadRhsFn,:N_Vector) shlib
@c Int32 CVodeQuadReInit (Ptr{:None},:N_Vector) shlib
@c Int32 CVodeQuadSStolerances (Ptr{:None},:realtype,:realtype) shlib
@c Int32 CVodeQuadSVtolerances (Ptr{:None},:realtype,:N_Vector) shlib
@c Int32 CVodeSensInit (Ptr{:None},:Int32,:Int32,:CVSensRhsFn,Ptr{:N_Vector}) shlib
@c Int32 CVodeSensInit1 (Ptr{:None},:Int32,:Int32,:CVSensRhs1Fn,Ptr{:N_Vector}) shlib
@c Int32 CVodeSensReInit (Ptr{:None},:Int32,Ptr{:N_Vector}) shlib
@c Int32 CVodeSensSStolerances (Ptr{:None},:realtype,Ptr{:realtype}) shlib
@c Int32 CVodeSensSVtolerances (Ptr{:None},:realtype,Ptr{:N_Vector}) shlib
@c Int32 CVodeSensEEtolerances (Ptr{:None},) shlib
@c Int32 CVodeQuadSensInit (Ptr{:None},:CVQuadSensRhsFn,Ptr{:N_Vector}) shlib
@c Int32 CVodeQuadSensReInit (Ptr{:None},Ptr{:N_Vector}) shlib
@c Int32 CVodeQuadSensSStolerances (Ptr{:None},:realtype,Ptr{:realtype}) shlib
@c Int32 CVodeQuadSensSVtolerances (Ptr{:None},:realtype,Ptr{:N_Vector}) shlib
@c Int32 CVodeQuadSensEEtolerances (Ptr{:None},) shlib
@c None CVodeQuadFree (Ptr{:None},) shlib
@c None CVodeSensFree (Ptr{:None},) shlib
@c None CVodeQuadSensFree (Ptr{:None},) shlib
@c Int32 CVodeSetQuadErrCon (Ptr{:None},:Int32) shlib
@c Int32 CVodeSetSensDQMethod (Ptr{:None},:Int32,:realtype) shlib
@c Int32 CVodeSetSensErrCon (Ptr{:None},:Int32) shlib
@c Int32 CVodeSetSensMaxNonlinIters (Ptr{:None},:Int32) shlib
@c Int32 CVodeSetSensParams (Ptr{:None},Ptr{:realtype},Ptr{:realtype},Ptr{:Int32}) shlib
@c Int32 CVodeSetQuadSensErrCon (Ptr{:None},:Int32) shlib
@c Int32 CVodeSensToggleOff (Ptr{:None},) shlib
@c Int32 CVodeGetQuad (Ptr{:None},Ptr{:realtype},:N_Vector) shlib
@c Int32 CVodeGetQuadDky (Ptr{:None},:realtype,:Int32,:N_Vector) shlib
@c Int32 CVodeGetSens (Ptr{:None},Ptr{:realtype},Ptr{:N_Vector}) shlib
@c Int32 CVodeGetSens1 (Ptr{:None},Ptr{:realtype},:Int32,:N_Vector) shlib
@c Int32 CVodeGetSensDky (Ptr{:None},:realtype,:Int32,Ptr{:N_Vector}) shlib
@c Int32 CVodeGetSensDky1 (Ptr{:None},:realtype,:Int32,:Int32,:N_Vector) shlib
@c Int32 CVodeGetQuadSens (Ptr{:None},Ptr{:realtype},Ptr{:N_Vector}) shlib
@c Int32 CVodeGetQuadSens1 (Ptr{:None},Ptr{:realtype},:Int32,:N_Vector) shlib
@c Int32 CVodeGetQuadSensDky (Ptr{:None},:realtype,:Int32,Ptr{:N_Vector}) shlib
@c Int32 CVodeGetQuadSensDky1 (Ptr{:None},:realtype,:Int32,:Int32,:N_Vector) shlib
@c Int32 CVodeGetQuadNumRhsEvals (Ptr{:None},Ptr{:Clong}) shlib
@c Int32 CVodeGetQuadNumErrTestFails (Ptr{:None},Ptr{:Clong}) shlib
@c Int32 CVodeGetQuadErrWeights (Ptr{:None},:N_Vector) shlib
@c Int32 CVodeGetQuadStats (Ptr{:None},Ptr{:Clong},Ptr{:Clong}) shlib
@c Int32 CVodeGetSensNumRhsEvals (Ptr{:None},Ptr{:Clong}) shlib
@c Int32 CVodeGetNumRhsEvalsSens (Ptr{:None},Ptr{:Clong}) shlib
@c Int32 CVodeGetSensNumErrTestFails (Ptr{:None},Ptr{:Clong}) shlib
@c Int32 CVodeGetSensNumLinSolvSetups (Ptr{:None},Ptr{:Clong}) shlib
@c Int32 CVodeGetSensErrWeights (Ptr{:None},Ptr{:N_Vector}) shlib
@c Int32 CVodeGetSensStats (Ptr{:None},Ptr{:Clong},Ptr{:Clong},Ptr{:Clong},Ptr{:Clong}) shlib
@c Int32 CVodeGetSensNumNonlinSolvIters (Ptr{:None},Ptr{:Clong}) shlib
@c Int32 CVodeGetSensNumNonlinSolvConvFails (Ptr{:None},Ptr{:Clong}) shlib
@c Int32 CVodeGetStgrSensNumNonlinSolvIters (Ptr{:None},Ptr{:Clong}) shlib
@c Int32 CVodeGetStgrSensNumNonlinSolvConvFails (Ptr{:None},Ptr{:Clong}) shlib
@c Int32 CVodeGetSensNonlinSolvStats (Ptr{:None},Ptr{:Clong},Ptr{:Clong}) shlib
@c Int32 CVodeGetQuadSensNumRhsEvals (Ptr{:None},Ptr{:Clong}) shlib
@c Int32 CVodeGetQuadSensNumErrTestFails (Ptr{:None},Ptr{:Clong}) shlib
@c Int32 CVodeGetQuadSensErrWeights (Ptr{:None},Ptr{:N_Vector}) shlib
@c Int32 CVodeGetQuadSensStats (Ptr{:None},Ptr{:Clong},Ptr{:Clong}) shlib
@c Int32 CVodeAdjInit (Ptr{:None},:Clong,:Int32) shlib
@c Int32 CVodeAdjReInit (Ptr{:None},) shlib
@c None CVodeAdjFree (Ptr{:None},) shlib
@c Int32 CVodeCreateB (Ptr{:None},:Int32,:Int32,Ptr{:Int32}) shlib
@c Int32 CVodeInitB (Ptr{:None},:Int32,:CVRhsFnB,:realtype,:N_Vector) shlib
@c Int32 CVodeInitBS (Ptr{:None},:Int32,:CVRhsFnBS,:realtype,:N_Vector) shlib
@c Int32 CVodeReInitB (Ptr{:None},:Int32,:realtype,:N_Vector) shlib
@c Int32 CVodeSStolerancesB (Ptr{:None},:Int32,:realtype,:realtype) shlib
@c Int32 CVodeSVtolerancesB (Ptr{:None},:Int32,:realtype,:N_Vector) shlib
@c Int32 CVodeQuadInitB (Ptr{:None},:Int32,:CVQuadRhsFnB,:N_Vector) shlib
@c Int32 CVodeQuadInitBS (Ptr{:None},:Int32,:CVQuadRhsFnBS,:N_Vector) shlib
@c Int32 CVodeQuadReInitB (Ptr{:None},:Int32,:N_Vector) shlib
@c Int32 CVodeQuadSStolerancesB (Ptr{:None},:Int32,:realtype,:realtype) shlib
@c Int32 CVodeQuadSVtolerancesB (Ptr{:None},:Int32,:realtype,:N_Vector) shlib
@c Int32 CVodeF (Ptr{:None},:realtype,:N_Vector,Ptr{:realtype},:Int32,Ptr{:Int32}) shlib
@c Int32 CVodeB (Ptr{:None},:realtype,:Int32) shlib
@c Int32 CVodeSetAdjNoSensi (Ptr{:None},) shlib
@c Int32 CVodeSetIterTypeB (Ptr{:None},:Int32,:Int32) shlib
@c Int32 CVodeSetUserDataB (Ptr{:None},:Int32,Ptr{:None}) shlib
@c Int32 CVodeSetMaxOrdB (Ptr{:None},:Int32,:Int32) shlib
@c Int32 CVodeSetMaxNumStepsB (Ptr{:None},:Int32,:Clong) shlib
@c Int32 CVodeSetStabLimDetB (Ptr{:None},:Int32,:Int32) shlib
@c Int32 CVodeSetInitStepB (Ptr{:None},:Int32,:realtype) shlib
@c Int32 CVodeSetMinStepB (Ptr{:None},:Int32,:realtype) shlib
@c Int32 CVodeSetMaxStepB (Ptr{:None},:Int32,:realtype) shlib
@c Int32 CVodeSetQuadErrConB (Ptr{:None},:Int32,:Int32) shlib
@c Int32 CVodeGetB (Ptr{:None},:Int32,Ptr{:realtype},:N_Vector) shlib
@c Int32 CVodeGetQuadB (Ptr{:None},:Int32,Ptr{:realtype},:N_Vector) shlib
@c Ptr{:None} CVodeGetAdjCVodeBmem (Ptr{:None},:Int32) shlib
@c Int32 CVodeGetAdjY (Ptr{:None},:realtype,:N_Vector) shlib
@ctypedef CVadjCheckPointRec Void
@c Int32 CVodeGetAdjCheckPointsInfo (Ptr{:None},Ptr{:CVadjCheckPointRec}) shlib
@c Int32 CVodeGetAdjDataPointHermite (Ptr{:None},:Int32,Ptr{:realtype},:N_Vector,:N_Vector) shlib
@c Int32 CVodeGetAdjDataPointPolynomial (Ptr{:None},:Int32,Ptr{:realtype},Ptr{:Int32},:N_Vector) shlib
@c Int32 CVodeGetAdjCurrentCheckPoint (Ptr{:None},Ptr{Ptr{:None}}) shlib

# header: /usr/local/include/cvodes/cvodes_impl.h
@ctypedef CVadjMem Ptr{:Void}
@ctypedef CkpntMem Ptr{:Void}
@ctypedef DtpntMem Ptr{:Void}
@ctypedef CVodeBMem Ptr{:Void}
@ctypedef cvaIMMallocFn Ptr{:Void}
@ctypedef cvaIMFreeFn Ptr{:Void}
@ctypedef cvaIMGetYFn Ptr{:Void}
@ctypedef cvaIMStorePntFn Ptr{:Void}
@ctypedef HermiteDataMem Ptr{:Void}
@ctypedef PolynomialDataMem Ptr{:Void}
@c Int32 cvEwtSet (:N_Vector,:N_Vector,Ptr{:None}) shlib
@c None cvProcessError (:CVodeMem,:Int32,Ptr{:Uint8},Ptr{:Uint8},Ptr{:Uint8}) shlib
@c None cvErrHandler (:Int32,Ptr{:Uint8},Ptr{:Uint8},Ptr{:Uint8},Ptr{:None}) shlib
@c Int32 cvSensRhsWrapper (:CVodeMem,:realtype,:N_Vector,:N_Vector,Ptr{:N_Vector},Ptr{:N_Vector},:N_Vector,:N_Vector) shlib
@c Int32 cvSensRhs1Wrapper (:CVodeMem,:realtype,:N_Vector,:N_Vector,:Int32,:N_Vector,:N_Vector,:N_Vector,:N_Vector) shlib
@c Int32 cvSensRhsInternalDQ (:Int32,:realtype,:N_Vector,:N_Vector,Ptr{:N_Vector},Ptr{:N_Vector},Ptr{:None},:N_Vector,:N_Vector) shlib
@c Int32 cvSensRhs1InternalDQ (:Int32,:realtype,:N_Vector,:N_Vector,:Int32,:N_Vector,:N_Vector,Ptr{:None},:N_Vector,:N_Vector) shlib

# header: /usr/local/include/cvodes/cvodes_spbcgs.h
@ctypedef CVSpilsPrecSetupFnB Ptr{:Void}
@ctypedef CVSpilsPrecSolveFnB Ptr{:Void}
@ctypedef CVSpilsJacTimesVecFnB Ptr{:Void}
@c Int32 CVSpilsSetPrecTypeB (Ptr{:None},:Int32,:Int32) shlib
@c Int32 CVSpilsSetGSTypeB (Ptr{:None},:Int32,:Int32) shlib
@c Int32 CVSpilsSetEpslinB (Ptr{:None},:Int32,:realtype) shlib
@c Int32 CVSpilsSetMaxlB (Ptr{:None},:Int32,:Int32) shlib
@c Int32 CVSpilsSetPreconditionerB (Ptr{:None},:Int32,:CVSpilsPrecSetupFnB,:CVSpilsPrecSolveFnB) shlib
@c Int32 CVSpilsSetJacTimesVecFnB (Ptr{:None},:Int32,:CVSpilsJacTimesVecFnB) shlib
@c Int32 CVSpbcgB (Ptr{:None},:Int32,:Int32,:Int32) shlib

# header: /usr/local/include/cvodes/cvodes_spgmr.h
@c Int32 CVSpgmrB (Ptr{:None},:Int32,:Int32,:Int32) shlib

# header: /usr/local/include/cvodes/cvodes_spils.h

# header: /usr/local/include/cvodes/cvodes_sptfqmr.h
@c Int32 CVSptfqmrB (Ptr{:None},:Int32,:Int32,:Int32) shlib


text:
mutable struct RRSIBPVAE <: AbstractIBPVAE
    prior
    N::Int
    k_max::Int
    x_dim::Int
    h_dim::Int
    encoder
    decoder
    rho
end

function Base.show(io::IO, vae::RRSIBPVAE)
    s = "RRSIBPVAE"
    s *= "\n  #params = $(numparams(vae))"
    print(io, s)
end

function (vae::RRSIBPVAE)(x, ::Val{true}; writer=nothing, epoch=nothing, iter=nothing, m::Int=10, τ=FT(0.1), kwargs...)
    # Roulette
    k_list = rand(value(vae.rho).lnpd, m)
    k_max = max(k_list...)
    if k_max > vae.k_max
        vae.k_max = k_max
    end

    # Compute weights for efficient multi-sample RRS
    rho = getρ(value(vae.rho), 1, k_max)
    one_minus_rho = 1 .- rho
    one_minus_rho_weighted = zeros(FT, k_max)
    # Multiple sample re-weighting
    for k in k_list
        one_minus_rho_weighted[1:k] .+= one_minus_rho[1:k]
    end
    one_minus_rho_weighted ./= m
    # NOTE: below takes use the independent property of KLs
    # .     cumsum due to the fact that R(i) includes all R(j<i)
    weight = AT(reverse(cumsum(reverse(one_minus_rho_weighted))))

    dist_nu, dist_gumbel, dist_A = vae.encoder(x, k_max)

    _, logpi = rand_logpi(dist_nu)
    logitpi = logit.(exp.(logpi))

    dist_Z = BatchGumbelBernoulliLogit(dist_gumbel.logitp .+ logitpi)
    logitZ = logitrand(dist_Z; τ=τ)

    L = zeros(FT, k_max)

    kl_A = if isdeep(vae)
        _kl_A = sum(kldiv(dist_A, BatchNormal(zero(FT), one(FT))); dims=2)
        L .+= vec(Array(value(_kl_A)))
        sum(_kl_A .* weight)
    else
        0
    end

    # Following Singh et al. (2017) use a multiplier of 1,000 for the KL term of nu to encourage adhering to the IBP prior.
    global_multiplier = size(x, 2) / vae.N * 1_000
    _kl_nu = sum(kldiv(dist_nu, BatchBeta(vae.prior.α, one(FT))); dims=2)
    _kl_nu = _kl_nu .* global_multiplier
    L .+= vec(Array(value(_kl_nu)))
    kl_nu = sum(_kl_nu .* weight)

    lq_Z = logpdflogit(dist_Z, logitZ; τ=τ)
    lp_Z = logpdflogit(BatchGumbelBernoulliLogit(logitpi), logitZ; τ=τ)
    _kl_Z = sum(lq_Z - lp_Z; dims=2)
    L .+= vec(Array(value(_kl_Z)))
    kl_Z = sum(_kl_Z .* weight)

    Z = sigm.(logitZ)
    code = isdeep(vae) ? rand(dist_A) .* Z : Z
    local lp_x
    for k = 1:k_max
        dist_x_k = vae.decoder(code[1:k,:])
        lp_x_k = sum(logpdf(dist_x_k, x))
        L[k] -= value(lp_x_k)
        if k == 1
            lp_x = lp_x_k * one_minus_rho_weighted[k]
        else
            lp_x += lp_x_k * one_minus_rho_weighted[k]
        end
    end

    L ./= size(x, 2)

    write_ibp_train_stats!(isdeep(vae), writer, iter, kl_A, kl_nu, kl_Z, lp_x, Z; k_list=k_list)

    # Update gradient for Rho
    dLdrho = MLToolkit.gradRR(i -> L[i], AutoGrad.value(vae.rho), k_list)
    updategrad!(AutoGrad.value(vae.rho), dLdrho)

    return (kl_nu + kl_Z + kl_A - lp_x) / size(x, 2)
end

function (vae::RRSIBPVAE)(x, ::Val{false}; m=10)
    k_list = rand(value(vae.rho).lnpd, 1_000)
    k = floor(Int, mean(k_list)) # mean
    return iwae_structured(vae, x, m, k)
end


text:
using NetworkModels
using Test

@testset "NetworkModels.jl" begin
    # Write your tests here.
end


text:
using DashBootstrapComponents, DashHtmlComponents

card_content = [
    dbc_cardheader("Card header"),
    dbc_cardbody([
        html_h5("Card title", className="card-title"),
        html_p("This is some card content that we'll reuse", className="card-text"),
    ]),
];

row_1 = dbc_row(
    [
        dbc_col(dbc_card(card_content, color="primary", outline=true)),
        dbc_col(dbc_card(card_content, color="secondary", outline=true)),
        dbc_col(dbc_card(card_content, color="info", outline=true)),
    ],
    className="mb-4",
);

row_2 = dbc_row(
    [
        dbc_col(dbc_card(card_content, color="success", outline=true)),
        dbc_col(dbc_card(card_content, color="warning", outline=true)),
        dbc_col(dbc_card(card_content, color="danger", outline=true)),
    ],
    className="mb-4",
);

row_3 = dbc_row([
    dbc_col(dbc_card(card_content, color="light", outline=true)),
    dbc_col(dbc_card(card_content, color="dark", outline=true)),
]);

cards = html_div([row_1, row_2, row_3]);

max_stars_repo_name: harisbal/dash-bootstrap-components

max_issues_repo_name: harisbal/dash-bootstrap-components

max_forks_repo_name: harisbal/dash-bootstrap-components


text:
# This file was generated by the Julia Swagger Code Generator
# Do not modify this file directly. Modify the swagger specification instead.

mutable struct IoK8sApiCoreV1SELinuxOptions <: SwaggerModel
    level::Any # spec type: Union{ Nothing, String } # spec name: level
    role::Any # spec type: Union{ Nothing, String } # spec name: role
    type::Any # spec type: Union{ Nothing, String } # spec name: type
    user::Any # spec type: Union{ Nothing, String } # spec name: user

    function IoK8sApiCoreV1SELinuxOptions(;level=nothing, role=nothing, type=nothing, user=nothing)
        o = new()
        validate_property(IoK8sApiCoreV1SELinuxOptions, Symbol("level"), level)
        setfield!(o, Symbol("level"), level)
        validate_property(IoK8sApiCoreV1SELinuxOptions, Symbol("role"), role)
        setfield!(o, Symbol("role"), role)
        validate_property(IoK8sApiCoreV1SELinuxOptions, Symbol("type"), type)
        setfield!(o, Symbol("type"), type)
        validate_property(IoK8sApiCoreV1SELinuxOptions, Symbol("user"), user)
        setfield!(o, Symbol("user"), user)
        o
    end
end # type IoK8sApiCoreV1SELinuxOptions

const _property_map_IoK8sApiCoreV1SELinuxOptions = Dict{Symbol,Symbol}(Symbol("level")=>Symbol("level"), Symbol("role")=>Symbol("role"), Symbol("type")=>Symbol("type"), Symbol("user")=>Symbol("user"))
const _property_types_IoK8sApiCoreV1SELinuxOptions = Dict{Symbol,String}(Symbol("level")=>"String", Symbol("role")=>"String", Symbol("type")=>"String", Symbol("user")=>"String")
Base.propertynames(::Type{ IoK8sApiCoreV1SELinuxOptions }) = collect(keys(_property_map_IoK8sApiCoreV1SELinuxOptions))
Swagger.property_type(::Type{ IoK8sApiCoreV1SELinuxOptions }, name::Symbol) = Union{Nothing,eval(Base.Meta.parse(_property_types_IoK8sApiCoreV1SELinuxOptions[name]))}
Swagger.field_name(::Type{ IoK8sApiCoreV1SELinuxOptions }, property_name::Symbol) =  _property_map_IoK8sApiCoreV1SELinuxOptions[property_name]

function check_required(o::IoK8sApiCoreV1SELinuxOptions)
    true
end

function validate_property(::Type{ IoK8sApiCoreV1SELinuxOptions }, name::Symbol, val)
end


text:
using BitFlags
using Test, Serialization

macro macrocall(ex)
    @assert Meta.isexpr(ex, :macrocall)
    ex.head = :call
    for i in 2:length(ex.args)
        ex.args[i] = QuoteNode(ex.args[i])
    end
    insert!(ex.args, 3, __module__)
    return esc(ex)
end

#@testset "Basic operation" begin
    # Inline definition
    @bitflag Flag1 flag1a flag1b flag1c
    @test Int(flag1a) == 1
    @test Int(flag1b) == 2
    @test Int(flag1c) == 4
    @test Flag1(1) == flag1a
    @test Flag1(2) == flag1b
    @test Flag1(4) == flag1c
    @test_throws ArgumentError Flag1(0)
    @test_throws ArgumentError Flag1(8)
    @test instances(Flag1) == (flag1a, flag1b, flag1c)

    # Block definition
    @bitflag Flag2 begin
        flag2a
        flag2b
        flag2c
    end
    @test Int(flag2a) == 1
    @test Int(flag2b) == 2
    @test Int(flag2c) == 4

    # Explicit numbering, inline
    @bitflag Flag3 flag3a=1 flag3b flag3c=8
    @test Int(flag3a) == 1
    @test Int(flag3b) == 2
    @test Int(flag3c) == 8
    @test typemin(Flag3) == flag3a
    @test typemax(Flag3) == flag3c
    @test Int(typemin(Flag3)) == 1
    @test Int(typemax(Flag3)) == 8
    @test_throws ArgumentError Flag3(4)

    # Zero-valued bit flag
    @bitflag Flag4 begin
        flag4a = 0
        flag4b
        flag4c
    end
    @test Int(flag4a) == 0
    @test Int(flag4b) == 1
    @test Int(flag4c) == 2
    @test Flag4(0) == flag4a

    # Mask operations
    @test Int(Flag1(7)) == 7
    @test Flag1(7) == flag1a | flag1b | flag1c
    @test Flag1(7) & flag1a == flag1a
    @test flag1a < flag1b < flag1c
    @test flag1a | flag1b < flag1c
#end

#@testset "Type properties" begin
    # Default integer typing
    @bitflag Flag5 flag5a flag5b
    @test typeof(Integer(flag5a)) == UInt32
    @test typeof(Flag5) == DataType
    @test typeof(flag5a) <: Flag5 <: BitFlag
    @test isbitstype(Flag5)
    @test isbits(flag5a)

    # Construct non-default
    @bitflag Flag6::UInt8 flag6a flag6b flag6c
    @test typeof(Integer(flag6a)) == UInt8
    @test UInt8(flag6a) === 0x01
    @test UInt16(flag6b) === 0x0002
    @test UInt128(flag6c) === 0x00000000000000000000000000000004

    # Explicit values of non-default types
    @bitflag Flag7::UInt8 flag7a=big"1" flag7b=UInt8(2)
    @test Integer(flag7a) === UInt8(1)
    @test Integer(flag7b) === UInt8(2)
#end

#@testset "Error conditions" begin

    @test_throws ArgumentError("no arguments given for BitFlag Foo"
                              ) @macrocall(@bitflag Foo)
    @test_throws ArgumentError("invalid base type for BitFlag Foo, "
                               * "::Float64; base type must be an unsigned "
                               * "integer primitive type"
                              ) @macrocall(@bitflag Foo::Float64 x=1.)

    # Require uniqueness
    @test_throws ArgumentError("values for BitFlag Foo are not unique"
                              ) @macrocall(@bitflag Foo x=1 y=1)
    @test_throws ArgumentError("values for BitFlag Foo are not unique"
                              ) @macrocall(@bitflag Foo x=0 y=0)

    # Explicit values must be powers of two
    @test_throws ArgumentError("invalid value for BitFlag Foo, _three = 3; "
                               * "values must be a positive power of 2"
                              ) @macrocall(@bitflag Foo _three = 3)

    # Values must be integers
    @test_throws ArgumentError("invalid value for BitFlag Foo, _zero = \"zero\"; "
                               * "values must be unsigned integers"
                              ) @macrocall(@bitflag Foo _zero="zero")
    @test_throws ArgumentError("invalid value for BitFlag Foo, _zero = '0'; "
                               * "values must be unsigned integers"
                              ) @macrocall(@bitflag Foo _zero='0')
    @test_throws ArgumentError("invalid value for BitFlag Foo, _zero = 0.5; "
                               * "values must be unsigned integers"
                              ) @macrocall(@bitflag Foo _zero=0.5)

    # Names must be valid identifiers
    @test_throws ArgumentError("""invalid argument for BitFlag Foo: if x
                                      1
                                  else
                                      2
                                  end"""
                              ) @macrocall(@bitflag Foo x ? 1 : 2)
    @test_throws ArgumentError("invalid argument for BitFlag Foo: 1 = 2"
                              ) @macrocall(@bitflag Foo 1=2)

    # Disallow value overflow
    @test_throws ArgumentError("overflow in value \"y\" of BitFlag Foo"
                          ) @macrocall(@bitflag Foo x=(UInt32(1) << 31) y)

#end

#@testset "Input/Output" begin
    @bitflag Flag9::UInt64 flag9a=1 flag9b=256
    @bitflag Flag10::UInt8 flag10a=128 flag10b=2
    let io = IOBuffer()
        write(io, flag10a | flag10b)
        write(io, flag9b)
        write(io, flag10a)
        seekstart(io)
        @test read(io, Flag10) == flag10a | flag10b
        @test read(io, Flag9) == flag9b
        @test read(io, Flag10) == flag10a
    end
    let io = IOBuffer()
        serialize(io, flag9a)
        seekstart(io)
        @test deserialize(io) === flag9a
    end
#end

#@testset "String representations" begin
    @bitflag FilePerms::UInt8 NONE=0 READ=4 WRITE=2 EXEC=1
    # In submodule and not imported to test prefix printing
    module SubModule
        using ..BitFlags
        @bitflag Bits::UInt8 BIT_ONE BIT_TWO BIT_FOUR BIT_EIGHT # No-zero case
    end

    @test string(FilePerms) == "FilePerms"
    @test string(SubModule.Bits) == "Main.SubModule.Bits"
    @test string(NONE) == "NONE"
    @test string(SubModule.BIT_ONE) == "BIT_ONE"
    @test repr("text/plain", FilePerms) ==
        """BitFlag $(string(FilePerms)):
           NONE = 0x00
           EXEC = 0x01
           WRITE = 0x02
           READ = 0x04"""
    @test repr("text/plain", SubModule.Bits) ==
        """BitFlag Main.SubModule.Bits:
           BIT_ONE = 0x01
           BIT_TWO = 0x02
           BIT_FOUR = 0x04
           BIT_EIGHT = 0x08"""
    @test repr(EXEC) == "EXEC::FilePerms = 0x01"
    @test repr(SubModule.BIT_ONE) == "BIT_ONE::Main.SubModule.Bits = 0x01"
    @test repr(EXEC | READ) == "(READ | EXEC)::FilePerms = 0x05"
    @test repr(SubModule.BIT_ONE | SubModule.BIT_EIGHT) == "(BIT_EIGHT | BIT_ONE)::Main.SubModule.Bits = 0x09"
    @test repr(NONE | READ) == "READ::FilePerms = 0x04"

    let io = IOBuffer(),
        ioc = IOContext(io, :compact => true, :module => Main),
        iof = IOContext(io, :compact => false, :module => Main)
        # Explicit :compact => false required for consistency across Julia versions

        stringf(x) = (show(iof, MIME"text/plain"(), x); String(take!(io)))
        stringc(x) = (show(ioc, MIME"text/plain"(), x); String(take!(io)))

        @test stringf(BitFlag) == "BitFlag"
        @test stringc(BitFlag) == "BitFlag"

        @test stringf(Union{FilePerms, SubModule.Bits}) == "Union{FilePerms, Main.SubModule.Bits}"
        @test stringc(Union{FilePerms, SubModule.Bits}) == "Union{FilePerms, Bits}"

        @test stringc(NONE) == "NONE"
        @test stringc(SubModule.BIT_ONE) == "BIT_ONE"
        @test stringc(EXEC | READ) == "READ|EXEC"
        @test stringc(SubModule.BIT_ONE | SubModule.BIT_EIGHT) == "BIT_EIGHT|BIT_ONE"

        # Handle exceptional cases where set bits are not in the allowed set
        @test stringf(reinterpret(SubModule.Bits, 0x00)) == "reinterpret(Main.SubModule.Bits, 0x00)::Main.SubModule.Bits = 0x00"
        @test stringc(reinterpret(SubModule.Bits, 0x00)) == "0x00"

        @test stringf(reinterpret(SubModule.Bits, 0x11)) == "(BIT_ONE | reinterpret(Main.SubModule.Bits, 0x10))::Main.SubModule.Bits = 0x11"
        @test stringc(reinterpret(SubModule.Bits, 0x11)) == "BIT_ONE|0x10"

        @test stringf(reinterpret(FilePerms, 0x08)) == "reinterpret(FilePerms, 0x08)::FilePerms = 0x08"
        @test stringc(reinterpret(FilePerms, 0x08)) == "0x08"

        @test stringf(reinterpret(FilePerms, 0xff)) == "(READ | WRITE | EXEC | reinterpret(FilePerms, 0xf8))::FilePerms = 0xff"
        @test stringc(reinterpret(FilePerms, 0xff)) == "READ|WRITE|EXEC|0xf8"
    end
#end


text:
using SciMLTutorials
target = ARGS[1]
if isdir(target)
    if !isfile(joinpath(target, "Project.toml"))
        error("Cannot weave folder $(target) without Project.toml!")
    end
    println("Weaving the $(target) folder")
    SciMLTutorials.weave_folder(target)
elseif isfile(target)
    folder = dirname(target)
    file = basename(target)
    println("Weaving $(folder)/$(file)")
    SciMLTutorials.weave_file(folder, file)
else
    error("Unable to find weaving target $(target)!")
end


text:
## Graph type definition ##

"""
    Graph{N,E,D,M}

Type for storing graphs. Its fields have mostly self-explanatory names
(field marked with an asterisk is for internal use only):

* `nodes`      — vector with nodes in the graph
* `nodeattrs`  — dictionary with node attributes
* `edgemaxid`* — largest edge id
* `edgecount`  — count of edges in the graph
* `edgeattrs`  — dictionary with edge attributes

There are four type parameters:

* `N` — integer type used for node ids
* `E` — integer type used for edge ids, `Nothing` can also be used to save space
* `D` — singleton type indicating whether the graph is directed
    - `Directed`   — directed graph (default)
    - `Undirected` — undirected graph
* `M` — singleton type indicating whether multiple edges and self-loops
are allowed
    - `Multi`  — multiple edges and self-loops allowed, i.e. graph
    is a multigraph (default)
    - `Simple` — only single edges between pais of different nodes allowed,
    i.e. graph is simple

There is a couple of predefined type aliases you can use for dispatch in your
functions: `DirectedGraph` or `UndirectedGraph` (w.r.t graph directedness), and
`MultiGraph` or `SimpleGraph` (w.r.t. being a multigraph).

"""
mutable struct Graph{N<:Integer,E<:Union{Integer,Nothing},D<:DirParam,M<:MultiParam}
    nodes::Vector{Node{N,E}}
    nodeattrs::AttributeDict
    edgemaxid::E
    edgecount::Int
    edgeattrs::AttributeDict
end

const LightGraph = Graph{N,Nothing,D,M} where {N,D,M}

const DirectedGraph   = Graph{N,E,<:Directed,M} where {N,E,M}
const UndirectedGraph = Graph{N,E,<:Undirected,M} where {N,E,M}

const MultiGraph  = Graph{N,E,D,Multi} where {N,E,D}
const SimpleGraph = Graph{N,E,D,Simple} where {N,E,D}

Graph{N,E,D,M}() where {N,E,D,M} = Graph{N,E,D,M}(
        Node{N,E}[], AttributeDict(),          # nodes
        defval(E), zero(Int), AttributeDict()  # edges
    )

Graph{N,E,D,M}(g::Graph) where {N,E,D,M} = Graph{N,E,D,M}(
        g.nodes, g.nodeattrs,                  # nodes
        g.edgemaxid, g.edgecount, g.edgeattrs  # edges
    )

dirview(g::Graph{N,E,D_,M}, D::Type{<:DirParam}) where {N,E,D_,M} = Graph{N,E,D,M}(g)

"""
    isdirected(g::Graph)

Check whether graph is directed.
"""
isdirected(::DirectedGraph)   = true
isdirected(::UndirectedGraph) = false

"""
    ismultigraph(g::Graph)

Check whether graph is permitted to contain mutliple edges and/or selfloops.
"""
ismultigraph(::MultiGraph)  = true
ismultigraph(::SimpleGraph) = false

"""
    directed(g::Graph)

A view on the graph that makes it directed. Made without copying.
"""
directed(g::UndirectedGraph) = dirview(g, Forward)
directed(g::DirectedGraph)   = g

"""
    undirected(g::Graph)

A view on the graph that makes it undirected. Made without copying.
"""
undirected(g::DirectedGraph)   = dirview(g, Both)
undirected(g::UndirectedGraph) = g

"""
    reverse(g::Graph)

A view on the graph where each edge is reversed.
"""
reverse(g::Graph{N,E,Forward}) where {N,E} = dirview(g, Reverse)
reverse(g::Graph{N,E,Reverse}) where {N,E} = dirview(g, Forward)
reverse(g::UndirectedGraph) = g
adjoint(g::Graph) = reverse(g)


text:
module FindMarker
    export find_marker!

    using HorizonSideRobots
    include("../8/horizonside.jl")
    include("functional_robot.jl")

    """
    find_marker!()
    перемещает Робота, находящегося где-то на неораниченном поле без внутренних перегородок, в клетку с маркером
    """
    function find_marker!(robot)
        robot = interface_protected_robot(robot)
        rectangular_borders = interface_rectangular_borders(robot)

        function move_ifnomarker!(side)::Bool
            if robot.ismarker()
                return false
            end
            rectangular_borders.move!(side) # - return true
        end

        trajectories = interface_trajectories(robot)
        trajectories.spiral!(move_ifnomarker!)
    end
end

# ------- Исполняемая часть файла
using HorizonSideRobots
using .FindMarker
robot = Robot("9/example_9.2.sit"; animate=true)
FindMarker.find_marker!(robot)


text:
# Note that this script can accept some limited command-line arguments, run
# `julia build_tarballs.jl --help` to see a usage message.
using BinaryBuilder

name = "adwaita_icon_theme"
version = v"3.33.92"

# Collection of sources required to build adwaita-icon-theme
sources = [
    ArchiveSource("https://gitlab.gnome.org/GNOME/adwaita-icon-theme/-/archive/$(version)/adwaita-icon-theme-$(version).tar.bz2",
                  "9e2078bf9e4d28f2a921fa88159733fe83a1fd37f8cbd768a5de3b83f44f0973"),
]

# Bash recipe for building across all platforms
script = raw"""
cd $WORKSPACE/srcdir/adwaita-icon-theme-*/

# Install native gtk+3.0 so that we get `gtk-encode-symbolic-svg`
apk add gtk+3.0 librsvg
./autogen.sh --prefix=${prefix} --build=${MACHTYPE} --host=${target}
./configure --prefix=${prefix} --build=${MACHTYPE} --host=${target}
make -j${nproc}
make install
"""

# These are the platforms we will build for by default, unless further
# platforms are passed in on the command line
platforms = [AnyPlatform()]

# The products that we will ensure are always built
products = Product[
    FileProduct("share/icons", :adwaita_icons_dir),
]

# Dependencies that must be installed before this package can be built
dependencies = [
    Dependency("hicolor_icon_theme_jll"),
]

# Build the tarballs.
build_tarballs(ARGS, name, version, sources, script, platforms, products, dependencies)

max_stars_repo_name: sharanry/Yggdrasil

max_issues_repo_name: sharanry/Yggdrasil

max_forks_repo_name: sharanry/Yggdrasil


text:
# ==============================================================================
# Display Max-Plus scalars and matrices

"""
    mpstyle

Memory for saving the style of display for neutral and absorbing elements. By
default the display of ScicosLab will be used.
"""
global mpstyle = 1; # ScicosLab style

"""
    mp_change_display(style::Int)

Change the style of behavior of functions `Base.show()`:
- `-Inf` are displayed either with `ε` (style 2 or 3) or `.` symbols (style 1).
- `0` are displayed either with `e `(style 3) or '0' symbols (style 1 or 2).
- else: `-Inf` and `0` are displayed in Julia default sytle (style 0).

If this function is not called, by default the ScicosLab style will be used
(style 1).
"""
function mp_change_display(style::Int)
    global mpstyle = min(max(style, 0), 4)
end

#
name(::Type{MP}) = "Max-Plus "
name(::Type{MI}) = "Min-Plus "

# Base function for display a Max-Plus number.
function mpshow(io::IO, x::Trop{T}) where T
    if (mpstyle == 0)
        show(io, x.v)
    elseif x == zero(Trop{T})
        (mpstyle == 1 || mpstyle == 2) ? (@printf io ".") : (@printf io "ε")
    elseif x == one(Trop{T})
        (mpstyle == 1 || mpstyle == 3) ? (@printf io "0") : (@printf io "e")
    elseif x.v == trunc(x.v)
        (@printf io "%d" x.v)
    else
        show(io, x.v)
    end
end

function mpshow(io::IO, A::ArrTrop)
    if (size(A,2) == 1)
        print(io, size(A,1), "-element Max-Plus vector:\n")
    else
        print(io, size(A,1), '×', size(A,2), " Max-Plus dense matrix:\n")
    end
    pretty_table(io, A, tf = tf_borderless, noheader = true)
end

function mpshow(io::IO, A::LinearAlgebra.Transpose{Trop, ArrTrop})
    print(io, size(A,1), '×', size(A,2), " Max-Plus transposed dense matrix:\n")
    pretty_table(io, A, tf = tf_borderless, noheader = true)
end

function mpshow(io::IO, V::LinearAlgebra.Transpose{Trop, VecTrop})
    print(io, size(A,1), "-element Max-Plus transposed vector:\n")
    pretty_table(io, V, tf = tf_borderless, noheader = true)
end

# Note: we force the sparse display like done in Julia 1.5 because since Julia
# 1.6 sparse matrices are displayed like dense matrices with dots for zeros.
# This sounds weird since displayed huge sparse matrices take the same space
# than dense matrix.
mpshow(io::IO, S::SpaTrop) = show(io, S)

"""
    LaTeX(io::IO, A::Array{MP})

Base function for convert a Max-Plus dense matrix to a LaTeX formula. Symbols of
neutral and absorbing elements depends on mp_change_display(style).

# Examples
```julia-repl
julia> LaTeX(stdout, MP([4 3; 7 -Inf]))
\\left[
\\begin{array}{*{20}c}
4 & 3 \\\\
7 & . \\\\
\\end{array}
\\right]
```
"""
# FIXME: not a column-major traversal
function LaTeX(io::IO, A::ArrTrop)
    (@printf io "\\left[\n\\begin{array}{*{20}c}\n")
    for i in 1:size(A,1)
        for j in 1:size(A,2)
            if A[i,j] == zero(MP)
                if mpstyle == 0
                    (@printf io "-\\infty")
                elseif mpstyle == 3 || mpstyle == 4
                    (@printf io "\\varepsilon")
                else
                    (@printf io ".")
                end
            elseif A[i,j] == one(MP)
                if mpstyle == 2 || mpstyle == 4
                    (@printf io "e")
                else
                    (@printf io "%d" A[i,j].v)
                end
            elseif A[i,j].v == trunc(A[i,j].v)
                (@printf io "%d" A[i,j].v)
            else
                show(io, A[i,j].v)
            end
            if j < size(A, 2)
                (@printf io " & ")
            end
        end
        (@printf io " \\\\\n")
    end
    (@printf io "\\end{array}\n\\right]\n")
end

# Base function for displaying a Max-Plus sparse matrix.
LaTeX(io::IO, S::SpaTrop) = LaTeX(io, full(S))

function LaTeX(A::ArrTrop)
    s = "\\left[\n\\begin{array}{*{20}c}\n"
    for i in 1:size(A,1)
        for j in 1:size(A,2)
            if A[i,j] == zero(MP)
                if mpstyle == 0
                    s = s * "-\\infty"
                elseif mpstyle == 3 || mpstyle == 4
                    s = s * "\\varepsilon"
                else
                    s = s * "."
                end
            elseif A[i,j] == one(MP)
                if mpstyle == 2 || mpstyle == 4
                    s = s * "e"
                else
                    s = s * string(Int64(A[i,j].v))
                end
            elseif A[i,j].v == trunc(A[i,j].v)
                s = s * string(Int64(A[i,j].v))
            else
                s = s * string(A[i,j].v)
            end
            if j < size(A, 2)
                s = s * " & "
            end
        end
        s = s * " \\\\\n"
    end
    return s * "\\end{array}\n\\right]\n"
end

"""
    show(io::IO, x::MP)

Display a Max-Plus number depending on the currently set style:

# Examples
```julia-repl
julia> mp_change_display(0); mpeye(2,2)
2×2 Max-Plus dense array:
  0.0  -Inf
 -Inf   0.0

julia> mp_change_display(1); mpeye(2,2)
2×2 Max-Plus dense array:
 0  .
 .  0

julia> mp_change_display(2); mpeye(2,2)
2×2 Max-Plus dense array:
 e  .
 .  e

julia> mp_change_display(3); mpeye(2,2)
2×2 Max-Plus dense array:
 0  ε
 ε  0

julia> mp_change_display(4); mpeye(2,2)
2×2 Max-Plus dense array:
 e  ε
 ε  e
```
"""
# Called by pretty_table() when REPL shows a MP matrix.
# julia> [MP(1) MP(2); MP(3) MP(4)]
Base.show(io::IO, x::Trop) = mpshow(io, x)

# Called by REPL when showing a MP scalar.
# julia> MP(1)
function Base.show(io::IO, ::MIME"text/plain", x::MP)
    print(io, name(MP))
    mpshow(io, x)
end

function Base.show(io::IO, ::MIME"text/plain", x::MI)
    print(io, name(MI))
    mpshow(io, x)
end

# Called by the REPL through the display() method. This function fixes
# misaliged columns made by the default show() Julia. We use the package
# PrettyTables.
# julia> [MP(1) MP(2); MP(3) MP(4)]
Base.show(io::IO, ::MIME"text/plain", A::ArrTrop) = mpshow(io, A)

# Convert a Max-Plus dense matrix to a LaTeX formula. Symbols of
# neutral and absorbing elements depends on mp_change_display(style).
Base.show(io::IO, ::MIME"text/latex", A::ArrTrop) = LaTeX(io, A)

# Convert a Max-Plus sparse matrix to a LaTeX formula. Symbols of
# neutral and absorbing elements depends on mp_change_display(style).
Base.show(io::IO, ::MIME"text/latex", x::SpaTrop) = LaTeX(io, A)

Base.show(io::IO, ::MIME"text/plain", A::LinearAlgebra.Transpose{Trop, Matrix{Trop}}) = mpshow(io, A)
Base.show(io::IO, ::MIME"text/plain", V::LinearAlgebra.Transpose{Trop, Vector{Trop}}) = mpshow(io, V)


text:
module SFuncs

using Base: reinit_stdio
using ..MultiInterface

using ...Scruff
using ..Utils
using ..Operators

import ..Operators: Support

Operators.@import_operators()

include("sfuncs/dist/dist.jl")
include("sfuncs/dist/cat.jl")
include("sfuncs/dist/constant.jl")
include("sfuncs/dist/flip.jl")
include("sfuncs/dist/normal.jl")

include("sfuncs/score/score.jl")
include("sfuncs/score/hardscore.jl")
include("sfuncs/score/softscore.jl")
include("sfuncs/score/logscore.jl")
include("sfuncs/score/functionalscore.jl")
include("sfuncs/score/normalscore.jl")
include("sfuncs/score/parzen.jl")

include("sfuncs/util/extend.jl")

include("sfuncs/conddist/conditional.jl")
include("sfuncs/conddist/det.jl")
include("sfuncs/conddist/invertible.jl")
include("sfuncs/conddist/table.jl")
include("sfuncs/conddist/discretecpt.jl")
include("sfuncs/conddist/lineargaussian.jl")
include("sfuncs/conddist/CLG.jl")
include("sfuncs/conddist/separable.jl")
include("sfuncs/conddist/switch.jl")

include("sfuncs/compound/generate.jl")
include("sfuncs/compound/apply.jl")
include("sfuncs/compound/chain.jl")
include("sfuncs/compound/mixture.jl")
include("sfuncs/compound/network.jl")
include("sfuncs/compound/serial.jl")
include("sfuncs/compound/expander.jl")

include("sfuncs/op_impls/bp_ops.jl")
include("sfuncs/op_impls/basic_ops.jl")

end


text:
module LoopThrottleTest

using LoopThrottle
using Test

function fast_while_loop_throttle(i0, i_increment, imax, rate, min_sleep_time)
    i = i0
    nloops = 0
    elapsed = @elapsed begin
        @throttle i while i < imax
            i += i_increment
            nloops += 1
        end max_rate=rate min_sleep_time=min_sleep_time
    end
    elapsed, nloops, i
end

function fast_while_loop_test(i0, i_increment, imax, rate, min_sleep_time)
    elapsed, nloops, i = fast_while_loop_throttle(i0, i_increment, imax, rate, min_sleep_time)
    @test nloops == div(imax - i0, i_increment)
    @test i == i0 + nloops * i_increment
    @test elapsed ≈ i_increment * (nloops - 1) / rate atol = 1e-2
end

@testset "fast while loop" begin
    # compile
    fast_while_loop_throttle(0, 2, 2000, Inf, 0.01)

    # test
    fast_while_loop_test(0, 2, 2000, 1000., 0.01)
    fast_while_loop_test(6, 3, 3000, 2000., 0.01)
    fast_while_loop_test(0, 1, 10, 5., 0.01)
end

@testset "fast for loop" begin
    n = 1000
    function f()
        result = 0
        @throttle i for i = 0 : n
            result += i
        end max_rate = 2 * n
        result
    end

    result = f()
    @test result == div(n * (n + 1), 2)
    elapsed = @elapsed f()
    @test elapsed ≈ 0.5 atol = 5e-2
end

@testset "fencepost" begin
    function f()
        @throttle i for i = 0 : 1
            nothing
        end
    end
    f()
    elapsed = @elapsed f()
    @test elapsed ≈ 1. atol = 1e-2
end

@testset "slow for loop" begin
    n = 10
    iteration_time = 0.2

    f = function ()
        @throttle i for i = 1 : n
            sleep(iteration_time)
        end max_rate = 10.
    end

    g = function ()
        for i = 1 : n
            sleep(iteration_time)
        end
    end

    f()
    g()
    @test @elapsed(f()) ≈ @elapsed(g()) atol = 1e-2
end

end


text:
@testset "add_Ga" begin
    Gb_h = 0.0347
    Ga_m = 0.3
    @test @inferred add_Ga(Gb_h, Ga_m) == NamedTuple()
    #
    # here cannot @inferred the names in the return
    # TODO
    Ga = add_Ga(Gb_h, Ga_m, :N20 => 2, :CH4 => 4)
    @test keys(Ga) == (:Ga_N20, :Ga_CH4)
    @test ≈(Ga.Ga_N20, 0.0164, rtol = 1/100)
    @test ≈(Ga.Ga_CH4, 0.0105, rtol = 1/100)
    #
    Gam = add_Ga(Gb_h, missing, :N20 => 2, :CH4 => 4)
    @test keys(Gam) == (:Ga_N20, :Ga_CH4)
    @test all(ismissing.(values(Gam)))
    #
    # DataFrame variant
    dfo = DataFrame(Gb_h=SA[Gb_h, missing, 0.055], Ga_m = SA[0.3,0.3,0.3])
    df = copy(dfo)
    df2 = @inferred add_Gb!(df)
    @test isequal(df2, dfo)
    df2 = @inferred add_Ga!(df, :N20 => 2, :CH4 => 4)
    @test df2 === df # mutating
    @test propertynames(df2)[end-1:end] == [:Ga_N20, :Ga_CH4]
    @test df2.Ga_N20[1] == Ga.Ga_N20
end

@testset "compute_Ram" begin
    zr, zh, d = thal.zr, thal.zh, 0.7*thal.zh
    ustar, wind = tha48[24, [:ustar, :wind]]
    df = copy(tha48)
    stability_correction!(df; z=zr, d)
    z0m = roughness_parameters(Roughness_wind_profile(), df; zh, zr, psi_m = df.psi_m).z0m
    Ram_r = @inferred compute_Ram(ResistanceWindZr(), ustar, wind)
    @test Ram_r ≈ 6.31 rtol = 1/100 # TODO check with R
    Ram_p =  @inferred compute_Ram(
        ResistanceWindProfile(), ustar; zr, d, z0m, psi_h = df.psi_h[24])
    @test Ram_p ≈ 5.41 rtol = 1/100 # TODO check with R
    #
    #Gb_h, Gb_CO2 = compute_Gb!(df, Thom1972())[24, [:Gb_h, :Gb_CO2]]
    @inferred compute_Ram!(df, ResistanceWindZr())
    #@test all(propertynames(df)[(end+1-length(Ram)):end] .== keys(Ram))
    @test propertynames(df)[end] == :Ra_m
    @test df.Ra_m[24] == Ram_r
    @inferred compute_Ram!(df, ResistanceWindProfile(); zr, d, z0m)
    @test df.Ra_m[24] == Ram_p
    Ram_skalar = df.Ra_m
    #
    # test zr as a vecttor
    df[!,:zri] .= zr
    df.zri[1] = zr/2
    @inferred compute_Ram!(df, ResistanceWindProfile(); zr=df.zri, d, z0m)
    @test df.Ra_m[2:end] == Ram_skalar[2:end]
    @test df.Ra_m[1] != Ram_skalar[1]
end

@testset "roughness_length_heat" begin
    z0m, kB_h = 1.2, 3.4
    z0h = @inferred roughness_length_heat(z0m, kB_h)
    @test z0h ≈ 0.0400 rtol = 1e-2
end

@testset "aerodynamic_conductance! only ustar and wind" begin
    # using only ustar and wind
    dfo = DataFrame(ustar = [0.1,missing,0.3], wind = [3.4, 2.8, 3.3])
    df = copy(dfo)
    @inferred aerodynamic_conductance!(df; Gb_model=Thom1972())
    # missing due to missing zr
    @test all(ismissing.(df.psi_h))
    @test all(ismissing.(df.psi_m))
    @test all(ismissing.(df.z0h))
    @inferred aerodynamic_conductance!(
        df; Gb_model=Thom1972(), zr = thal.zr, zh = thal.zh,
        stab_formulation = NoStabilityCorrection())
    @test all(iszero.(df.psi_h))
    @test all(iszero.(df.psi_m))
    @test all(ismissing.(df.z0h))
    @test propertynames(df)[(end-9):end] ==
        SA[:Gb_h, :Rb_h, :kB_h, :Gb_CO2, :Ra_m, :Ga_m, :Ga_h, :Ra_h, :Ga_CO2, :z0h]
    # compare with R
    @test all(isapproxm.(df.Gb_h, (0.0347, missing, 0.0723), rtol = 1e-3))
    @test all(isapproxm.(df.Ga_m, (0.00294, missing, 0.0273), rtol = 1e-3))
    @test all(isapproxm.(df.Ga_h, (0.00271, missing, 0.0198), rtol = 1e-3))
end

@testset "aerodynamic_conductance! Gb_Thom" begin
    df = tha48[:,Not(Cols(:Gb_h, :Ga_h))]
    # test Ga_m by wind_profile
    @inferred aerodynamic_conductance!(
        df; Gb_model=Thom1972(), Ram_model=ResistanceWindProfile(),
        zh = thal.zh, zr = thal.zr)
    @test propertynames(df)[(end-9):end] ==
        SA[:Gb_h, :Rb_h, :kB_h, :Gb_CO2, :Ra_m, :Ga_m, :Ga_h, :Ra_h, :Ga_CO2, :z0h]
end

@testset "aerodynamic_conductance! constant_kB1" begin
    # specify kB1 as a vecttor
    df = tha48[:,Not(Cols(:Gb_h, :Ga_h))]
    df[!,:kB_hi] .= 1.18
    df.kB_hi[1] = 1.18/2
    @inferred aerodynamic_conductance!(df; Gb_model=ConstantKB1(), kB_h = df.kB_hi,
        zh = thal.zh, zr = thal.zr)
    @test propertynames(df)[(end-9):end] ==
        SA[:Gb_h, :Rb_h, :kB_h, :Gb_CO2, :Ra_m, :Ga_m, :Ga_h, :Ra_h, :Ga_CO2, :z0h]
    @test all(isapproxm.(df.Gb_h[1:3], (0.375, 0.17, 0.167), rtol = 1e-2))
end

@testset "aerodynamic_conductance! Gb_Choudhury" begin
    leafwidth=0.1
    df = tha48[:,Not(Cols(:Gb_h, :Ga_h))]
    # test Ga_m by wind_profile
    @inferred aerodynamic_conductance!(
        df; Gb_model=Choudhury1988(), Ram_model=ResistanceWindProfile(),
        leafwidth, LAI=thal.LAI, zh = thal.zh, zr = thal.zr)
    @test propertynames(df)[(end-9):end] ==
        SA[:Gb_h, :Rb_h, :kB_h, :Gb_CO2, :Ra_m, :Ga_m, :Ga_h, :Ra_h, :Ga_CO2, :z0h]
    # compare with R
    @test all(isapproxm.(df.Gb_h[1:3], (0.157, 0.15, 0.151), rtol = 1e-2))
    @test all(isapproxm.(df.Ga_m[1:3], (0.0709, 0.0649, 0.0603), rtol = 1e-2))
end

@testset "aerodynamic_conductance! Gb_Su" begin
    Dl=0.01
    df = tha48[:,Not(Cols(:Gb_h, :Ga_h))]
    df[!,:LAI] .= thal.LAI
    df.LAI[1] = thal.LAI/2
    # test providing changing LAI
    @inferred aerodynamic_conductance!(df; Gb_model=Su2001(),
        Dl, LAI=df.LAI, zh=thal.zh, zr=thal.zr);
    @test propertynames(df)[(end-9):end] ==
        SA[:Gb_h, :Rb_h, :kB_h, :Gb_CO2, :Ra_m, :Ga_m, :Ga_h, :Ra_h, :Ga_CO2, :z0h]
    # compare with R
    @test all(isapproxm.(df.Gb_h[1:3], (0.202, 0.177, 0.167), rtol = 1e-2))
    @test all(isapproxm.(df.Ga_m[1:3], (0.0693, 0.0538, 0.0507), rtol = 1e-2))
end


text:
using LazyTaylorSeries
using Test

@testset "Basic usage" begin
    t  = Taylor1((t, i) -> (i == 1) ? 1.0 : 0.0, Float64[])  # define variable
    t2 = Taylor1((t, i) -> (i == 1) ? 1.0 : 0.0, Dict{Int,Float64}())  # define variable

    @test t[0] == 0 == t2[0]
    @test t[1] == 1 == t2[1]
    @test t[2] == 0 == t2[2]
    @test t[100] == 0 == t2[100]

    s = 1 + t
    @test (s[0], s[1], s[2], s[3]) == (1, 1, 0, 0)

    s = t + t
    @test (s[0], s[1], s[2], s[3]) == (0, 2, 0, 0)

    s = t * t
    @test (s[0], s[1], s[2], s[3]) == (0, 0, 1, 0)

    s = t^2
    @test (s[0], s[1], s[2], s[3]) == (0, 0, 1, 0)

    s = (1 + t)^2
    @test (s[0], s[1], s[2], s[3]) == (1, 2, 1, 0)

end

max_stars_repo_name: dpsanders/LazyTaylorSeries

max_issues_repo_name: dpsanders/LazyTaylorSeries

max_forks_repo_name: dpsanders/LazyTaylorSeries


text:
"""
lwca.jl

Linearized WCA (LWCA) reference system

References:

* A. Meyer et al. CHemical Physics 49 (1980) 147-152
"""

immutable LWCA{T <: IndependentReferenceSystem} <: AbstractWCA
  trial::T
  temp::Float64
  optim::Symbol
  struct::Symbol
end

function LWCA(trial::IndependentReferenceSystem, temp::Real; optim::Symbol = :default, struct::Symbol = :linear)
  @assert optim in [:default, :Sandler] "LWCA: unsupported method of optimization"
  @assert struct in [:linear, :full] "LWCA: unsupported method of structure calculation"
  LWCA(trial, convert(Float64, temp), optim, struct)
end

immutable OptimizedLWCA{T <: IndependentReferenceSystem} <: AbstractOptimizedWCA
  trial::T
  temp::Float64
  rmin::Matrix{Float64}
  C₁::Matrix{Float64}
  C₂::Matrix{Float64}
  K₁::Matrix{Float64}
  K₂::Matrix{Float64}
  residue::Float64 # residue in optimization
end

function TPTSystem(lwca::LWCA{AHS}, pert::Perturbation; kwargs...)
  N::Int = ncomp(lwca)

  T::Float64 = temperature(lwca)
  β::Float64 = 1 / (kB * T)

  approx::Symbol = lwca.trial.approx
  c::Vector{Float64} = composition(lwca)
  ρ₀::Vector{Float64} = numberdensity(lwca.trial)

  u::Array{Function,2} = pairpotential(pert)
  u′::Array{Function,2} = pairpotential_derivative(pert)

  rcore::Array{Float64,2} = coreradius(pert)
  rcut::Array{Float64,2} = cutoffradius(pert)
  rmin::Array{Float64,2} = pairpotential_minimizer(pert)
  umin::Array{Float64,2} = [ u[i,j](rmin[i,j]) for i in 1:N, j in 1:N ]

  # Repulsive-part of pairpotential
  u₀::Array{Function,2} = [ r -> ( r < rmin[i,j] ? u[i,j](r) - umin[i,j] : 0. )
                            for i in 1:N, j in 1:N ]

  function fopt(σ::Vector{Float64})::Float64
    ahs = AHS(approx, σ, ρ₀)
    σ₀ = hsdiameter(ahs)
    g₀ = contactvalue(ahs)
    g′₀ = contactgradient(ahs)

    residue::Float64 = 0

    for i in 1:N, j in 1:N
      i ≠ j && lwca.optim ≠ :Sandler && continue
      i > j && continue

      Y = g′₀[i,j] / g₀[i,j] * σ₀[i,j]

      X = (-2β * σ₀[i,j] * u′[i,j](σ₀[i,j]) + Y + 2) /
          (-β * σ₀[i,j] * u′[i,j](σ₀[i,j]) + Y + 2)

      if X ≤ 0
        residue += Inf
      else
        residue += c[i]*c[j] * abs(β * u₀[i,j](σ₀[i,j]) - log(X))
      end
    end

    return residue
  end

  Δσ = lwca.optim == :Sander ? 0.05 : 0.01

  opt = NLopt.Opt(:LN_BOBYQA, N)
  NLopt.min_objective!(opt, (σ, g) -> fopt(σ))
  NLopt.xtol_rel!(opt, 1e-5)
  NLopt.maxeval!(opt, 100)
  NLopt.lower_bounds!(opt, [ rcore[i,i] for i in 1:N ])
  NLopt.upper_bounds!(opt, [ 0.99rmin[i,i] for i in 1:N ])
  NLopt.initial_step!(opt, [ Δσ*rmin[i,i] for i in 1:N ])

  σ_init = [ hsdiameter_estimate(pert, u, rmin, T)[i,i] for i in 1:N ]
  (fmin, σ_wca, res) = NLopt.optimize(opt, σ_init)

  ahs = AHS(approx, σ_wca, ρ₀)

  σ₀ = hsdiameter(ahs)
  g₀ = contactvalue(ahs)
  g′₀ = contactgradient(ahs)
  Y = [ g′₀[i,j] / g₀[i,j] * σ₀[i,j] for i in 1:N, j in 1:N ]

  C₁ = Array{Float64,2}(N,N) # C-
  C₂ = Array{Float64,2}(N,N) # C+
  K₁ = Array{Float64,2}(N,N) # K-
  K₂ = Array{Float64,2}(N,N) # K+

  for i in 1:N, j in 1:N
    i > j && continue

    σᵢⱼ::Float64 = (σ_wca[i] + σ_wca[j]) / 2
    u₀_σ::Float64 = u₀[i,j](σᵢⱼ)
    u′_σ::Float64 = u′[i,j](σᵢⱼ)

    C₁[i,j] = C₁[j,i] = σᵢⱼ^2 * g₀[i,j] * exp(-β*u₀_σ)
    C₂[i,j] = C₂[j,i] = σᵢⱼ^2 * g₀[i,j] * (exp(-β*u₀_σ) - 1)

    K₁[i,j] = K₁[j,i] = -β*σᵢⱼ*u′_σ + Y[i,j] + 2
    K₂[i,j] = K₂[j,i] =
      -β * σᵢⱼ * u′_σ * exp(-β*u₀_σ) / (exp(-β*u₀_σ) - 1) + Y[i,j] + 2
  end

  if lwca.struct == :linear
    optwca = OptimizedLWCA(ahs, T, rmin, C₁, C₂, K₁, K₂, fmin)
  elseif lwca.struct == :full
    optwca = OptimizedWCA(ahs, T, rmin, u₀, fmin)
  end

  TPTSystem(optwca, pert; kwargs...)
end

function blipfunction(lwca::OptimizedLWCA)::Array{Function,2}
  @attach(lwca, C₁, C₂, K₁, K₂)

  N = ncomp(lwca.trial)

  T = temperature(lwca)
  β = 1 / (kB * T)

  σ::Array{Float64,2} = hsdiameter(lwca)

  ret = Array{Function,2}(N,N)

  for i in 1:N, j in 1:N
    i > j && continue

    C₁′ = K₁[i,j] * C₁[i,j] / σ[i,j]
    C₂′ = K₂[i,j] * C₂[i,j] / σ[i,j]

    R₁ = σ[i,j] / K₁[i,j]
    R₂ = - σ[i,j] / K₂[i,j]

    function C(r)::Float64
      if r < σ[i,j] - R₁
        0
      elseif r < σ[i,j]
        C₁[i,j] + C₁′ * (r - σ[i,j])
      elseif r < σ[i,j] + R₂
        C₂[i,j] + C₂′ * (r - σ[i,j])
      else
        0
      end
    end

    B(r)::Float64 = C(r) / r^2

    ret[i,j] = ret[j,i] = B
  end

  return ret
end

function paircorrelation(lwca::OptimizedLWCA)::Array{Function,2}
  N::Int = ncomp(lwca)

  g_hs::Array{Function,2} = paircorrelation(lwca.trial)
  B::Array{Function,2} = blipfunction(lwca)

  ret = Array{Function,2}(N,N)

  for i in 1:N, j in 1:N
    i > j && continue

    g(r)::Float64 = g_hs[i,j](r) + B[i,j](r)

    ret[i,j] = ret[j,i] = g
  end

  return ret
end

function structurefactor(lwca::OptimizedLWCA)::Array{Function,2}
  @attach(lwca, C₁, C₂, K₁, K₂)

  N::Int = ncomp(lwca)
  c::Vector{Float64} = composition(lwca)

  S₀::Array{Function,2} = structurefactor(lwca.trial)
  ρ::Float64 = totalnumberdensity(lwca)

  ret = Array{Function,2}(N,N)

  for i in 1:N, j in 1:N
    i > j && continue

    σ::Float64 = hsdiameter(lwca)[i,j]

    B̃(q)::Float64 = -4π / factorial(3) * q*σ * sin(q*σ) / (q*σ) *
                    ( σ*C₁[i,j] / K₁[i,j]^2 - σ*C₂[i,j] / K₂[i,j]^2 )

    # S(q) = S₀[i,j](q) + ρ * √(c[i]*c[j]) * B̃(q)
    S(q) = S₀[i,j](q) / (1 - ρ * √(c[i]*c[j]) * S₀[i,j](q) * B̃(q))

    ret[i,j] = ret[j,i] = S
  end

  return ret
end


text:
type LabelerFF
    word_f
    tag_f
    label_f
    nonlinear
    W
end

const none = 46

function setheads!{T}(s::State{T})
    @assert isfinal(s)
    s.head = nullstate(T)
    stack = [s]
    while !isempty(stack)
        s = pop!(stack)
        if !isnull(s.lchild)
            s.lchild.head = s
            push!(stack, s.lchild)
        end
        st = s
        while !isnull(st.lsibl.lchild)
            st.lsibl.lchild.head = s
            push!(stack, st.lsibl.lchild)
            st = st.lsibl
        end
        if !isnull(s.rchild)
            s.rchild.head = s
            push!(stack, s.rchild)
        end
        st = s
        while !isnull(st.rsibl.rchild)
            st.rsibl.rchild.head = s
            push!(stack, st.rsibl.rchild)
            st = st.rsibl
        end
    end
end

type LExample
    w::Vector{Int}
    t::Union{Matrix{Float32},Vector{Int}}
    l::Vector{Int}
    target::Int
end

function gen_trainexample(s::State)
    sent = s.tokens
    labelsize = targetsize(s.parser.labeler)
    w, t = genfeatures(s)
    l = Int[] # [label(s0 -> s0head), label(s0head -> s0head's head)]
    push!(l, sent[s.top].label)
    push!(l, tokenat(s, s.head).label)
    push!(l, tokenat(s, s.head.rchild).label)
    push!(l, tokenat(s, s.head.lchild).label)
    push!(l, tokenat(s, s.head.head.rchild).label)
    push!(l, tokenat(s, s.head.head.lchild).label)
    target = sent[s.top].label
    LExample(w, t, l, target)
end

targetsize(m::LabelerFF) = size(m.W[end].w.data)[1]

function initmodel!(parser::DepParser, labeler::Type{LabelerFF}, embed,
    topktags, nonlinear, sparsesizes, embedsizes, hiddensizes; inherit=false)
    T = Float32
    if inherit
        word_f = Embedding(copy(parser.model.word_f.ws), IntSet())
        tag_f  = Embedding(copy(parser.model.tag_f.ws), IntSet())
        label_f = Embedding(copy(parser.model.label_f.ws), IntSet())
    else
        if embed == ""
            info("USING EMBEDDINGS WITH UNIFORM DISTRIBUTION [-0.01, 0.01]")
            word_f = Embedding(T, length(parser.words), embedsizes[1])
        else
            info("USING EMBEDDINGS LOADED FROM $(embed)")
            word_f = Embedding(embed, Float32)
            scale!(word_f)
        end
        if topktags
            tag_f = Linear(T, length(parser.tags), embedsizes[2])
        else
            tag_f = Embedding(T, length(parser.tags), embedsizes[2])
        end
        label_f = Embedding(T, length(parser.labels), embedsizes[2])
    end
    indim = sum(sparsesizes .* embedsizes)
    outdim = length(parser.labels) - 1 # exclude NONE
    W = [myLinear(T, indim, hiddensizes[1]),
         myLinear(T, hiddensizes[1], outdim)]
    parser.labeler = LabelerFF(word_f, tag_f, label_f, nonlinear, W)
    info("INPUT: [S^word,S^tag,S^label] = ", sparsesizes)
    info("EMBED DIMS: [word,tag,label] = ", embedsizes)
    info("HIDDEN LAYER: ", hiddensizes)
    info("OUTPUT DIM: ", outdim)
    info("NONLINEAR: ", nonlinear)
end

# TODO: make State have id field
# to tell where the State is in a batch
# called from expand(::State ::Int)
@compat function (m::LabelerFF){T}(s::State{T}, act::Int)
    Var([0f0])
end

@compat function (m::LabelerFF)(batch::AbstractVector, istrain=true)
    wordvec, tagvec, labelvec = [], [], []
    for s in batch
        push!(wordvec, s[1])
        push!(tagvec, s[2])
        push!(labelvec, s[3])
    end
    wordmat = m.word_f(Var(hcat(wordvec...)))
    tagmat = m.tag_f(Var(hcat(tagvec...)))
    labelmat = m.label_f(Var(hcat(labelvec...)))
    x = concat(1, wordmat, tagmat, labelmat)
    x = m.nonlinear(m.W[1](x))
    x = dropout(x, 0.5, istrain)
    x = m.W[end](x)
    x
end

function genfeatures(s::State)
    # word, tag
    s0     = tokenat(s, s.top)
    s0l    = tokenat(s, s.lchild)
    s0ll1  = tokenat(s, s.lchild.lchild)
    s0ll2  = tokenat(s, s.lchild.lsibl.lchild)
    s0lr1  = tokenat(s, s.lchild.rchild)
    s0lr2  = tokenat(s, s.lchild.rsibl.rchild)
    s0l2   = tokenat(s, s.lsibl.lchild)
    s0l2l1 = tokenat(s, s.lsibl.lchild.lchild)
    s0l2l2 = tokenat(s, s.lsibl.lchild.lsibl.lchild)
    s0l2r1 = tokenat(s, s.lsibl.lchild.rchild)
    s0l2r2 = tokenat(s, s.lsibl.lchild.rsibl.rchild)
    s0l3   = tokenat(s, s.lsibl.lsibl.lchild)
    s0r    = tokenat(s, s.rchild)
    s0rr1  = tokenat(s, s.rchild.rchild)
    s0rr2  = tokenat(s, s.rchild.rsibl.rchild)
    s0rl1  = tokenat(s, s.rchild.lchild)
    s0rl2  = tokenat(s, s.rchild.lsibl.lchild)
    s0r2   = tokenat(s, s.rsibl.rchild)
    s0r2r1 = tokenat(s, s.rsibl.rchild.rchild)
    s0r2r2 = tokenat(s, s.rsibl.rchild.rsibl.rchild)
    s0r2l1 = tokenat(s, s.rsibl.rchild.lchild)
    s0r2l2 = tokenat(s, s.rsibl.rchild.lsibl.lchild)
    s0r3   = tokenat(s, s.rsibl.rsibl.rchild)

    s0h   = tokenat(s, s.head)
    s0h2  = tokenat(s, s.head.head)
    s0hl  = tokenat(s, s.head.lchild)
    s0hr  = tokenat(s, s.head.rchild)
    s0h2l  = tokenat(s, s.head.head.lchild)
    s0h2r  = tokenat(s, s.head.head.rchild)

    words = [s0.word, s0l.word, s0ll1.word, s0ll2.word, s0lr1.word, s0lr2.word,
             s0l2.word, s0l2l1.word, s0l2l2.word, s0l2r1.word, s0l2r2.word, s0l3.word,
             s0.word, s0r.word, s0rr1.word, s0rr2.word, s0rl1.word, s0rl2.word,
             s0r2.word, s0r2r1.word, s0r2r2.word, s0r2l1.word, s0r2l2.word, s0r3.word,
             s0h.word, s0h2.word, s0hl.word, s0hr.word, s0h2l.word, s0h2r.word]

    tags = [s0.tag, s0l.tag, s0ll1.tag, s0ll2.tag, s0lr1.tag, s0lr2.tag,
             s0l2.tag, s0l2l1.tag, s0l2l2.tag, s0l2r1.tag, s0l2r2.tag, s0l3.tag,
             s0.tag, s0r.tag, s0rr1.tag, s0rr2.tag, s0rl1.tag, s0rl2.tag,
             s0r2.tag, s0r2r1.tag, s0r2r2.tag, s0r2l1.tag, s0r2l2.tag, s0r3.tag,
             s0h.tag, s0h2.tag, s0hl.tag, s0hr.tag, s0h2l.tag, s0h2r.tag]

    # words = [s0.word, s0l.word, s0l2.word, s0l3.word, s02l.word, s02l2.word, s0r.word,
    #  s0r2.word, s0r3.word, s02r.word, s02r2.word, s0h.word, s0h2.word]
    # tags = [s0.tag, s0l.tag, s0l2.tag, s0l3.tag, s02l.tag, s02l2.tag, s0r.tag,
    #         s0r2.tag, s0r3.tag, s02r.tag, s02r2.tag, s0h.tag, s0h2.tag,]
    return words, tags
end

# function getchildren(heads::Array{Int})
#     res = map(_ -> [], heads)
#     for (c, h) in enumerate(heads)
#         push!(res[h], c)
#     end
#     res
# end

# function head2state{T}(s::State{T})
#     @assert isfinal(s) && s.top == 0
#     res = Array(State{T}, length(s.tokens))
#     while !isnull(s.prev)
#         s = s.prev
#         if s.top != 0 && !isdefined(res, s.top)
#             res[s.top] = s
#         end
#     end
#     res
# end

function assign!{T}(labeler::LabelerFF, ss::Vector{State{T}})
    # predicted labels
    blabels = map(s ->fill(none, length(s.tokens)), ss)
    bheads  = map(s -> heads(s), ss)
    stack = map(i -> (ss[i].prev ,i), 1:length(ss)) # (batches i'th tree, batchid)
    while !isempty(stack)
        bsize = min(length(stack), 10000)
        batch, stack = stack[1:bsize], stack[bsize+1:end]
        bfeat = map(batch) do tup
            s, bid = tup
            w, t = genfeatures(s)
    push!(l, tokenat(s, s.head.rchild).label)
    push!(l, tokenat(s, s.head.lchild).label)
    push!(l, tokenat(s, s.head.head.rchild).label)
    push!(l, tokenat(s, s.head.head.lchild).label)
            l = [blabels[bid][s.top],
                get(blabels[bid], bheads[bid][s.top], none),
                s.head.rchild.top != 0 ? blabels[bid][bheads[bid][s.head.rchild.top]] : none,
                s.head.lchild.top != 0 ? blabels[bid][bheads[bid][s.head.lchild.top]] : none,
                s.head.head.rchild.top != 0 ? blabels[bid][bheads[bid][s.head.head.rchild.top]] : none,
                s.head.head.lchild.top != 0 ? blabels[bid][bheads[bid][s.head.head.lchild.top]] : none]
            (w, t, l)
        end
        preds = argmax(labeler(bfeat, false).data, 1)
        for (i, (s, bid)) in enumerate(batch)
            blabels[bid][s.top] = preds[i]
            !isnull(s.lchild) && push!(stack, (s.lchild, bid))
            st = s
            while !isnull(st.lsibl.lchild)
                push!(stack, (st.lsibl.lchild, bid))
                st = st.lsibl
            end
            !isnull(s.rchild) && push!(stack, (s.rchild, bid))
            st = s
            while !isnull(st.rsibl.rchild)
                push!(stack, (st.rsibl.rchild, bid))
                st = st.rsibl
            end
        end
    end
    @assert all(la -> all(v -> v != none, la), blabels)
    blabels
end
# if  && blabels[bid][s.(arc).top] == none

# function assign!{T}(labeler::LabelerFF, ss::Vector{State{T}})
#     @assert all(isfinal, ss)
#     # predicted labels
#     blabels   = map(s -> fill(none, length(s.tokens)), ss)
#     bstates   = map(s -> head2state(s), ss) # no State with s.top == 0
#     bheads    = map(s -> heads(s), ss)
#     bchildren = map(heads -> getchildren(heads), bheads)
#     stack = map(i -> (1,i), 1:length(ss)) # (headid, batchid)
#     while !isempty(stack)
#         bsize = min(length(stack), length(ss))
#         batch, stack = stack[1:bsize], stack[bsize:end]
#         bfeat = map(batch) do tup
#             hid, bid = tup
#             s = bstates[bid][hid]
#             w, t = genfeatures(s)
#             l = [labels[bid][s.top],
#                  heads[bid][s.top] != 0 ? labels[bid][heads[bid][s.top]] : 1]
#             (w, t, l)
#         end
#         preds = argmax(parser.labeler(bfeat, false).data, 1)
#         for (i, (hid, bid)) in enumarate(batch)
#             blabels[bid][hid] = preds[i]
#             for childid in bchildren
#                 push!(stack, (childid, bid))
#             end
#         end
#     end
#     blabels
# end
#
typealias Doc Vector{Vector{Token}}
function train!{T}(::Type{LabelerFF}, parser::DepParser{T}, trainsents::Doc,
    testsents::Doc=Vector{Token}[]; embed="", batchsize=32, iter=20, progbar=true,
    nonlinear=tanh, sparsesizes=[30,30,6] ,embedsizes=[50,50,50], hiddensizes=[200],
    topktags=false, opt=MyAdaGrad(0.01), evaliter=100, outfile="parser.dat")
    info("WILL RUN $iter ITERATIONS")

    if !isdefined(parser, :model)
        initmodel!(parser, FeedForward, embed, topktags,
               nonlinear, sparsesizes, embedsizes, hiddensizes)
    end
    saver = ModelSaver(outfile)
    initmodel!(parser, LabelerFF, embed, topktags,
               nonlinear, sparsesizes, embedsizes, hiddensizes)

    projectivize!(trainsents, false)

    trainsamples = LExample[]
    oracletrees = State{T}[]
    for sent in trainsents
        s = State(sent, parser)
        exs = Array(State{T}, length(sent))
        # exs[i] should be the final State object
        # at which i'th token appeared at s0.
        while !isfinal(s)
            s = expandgold(s)[1]
            if s.top != 0
                exs[s.top] = s
            end
        end
        setheads!(s)
        push!(oracletrees, s)
        @assert all(i -> isdefined(exs, i), 1:length(exs))
        append!(trainsamples, map(gen_trainexample, exs))
    end
    @assert all(isfinal, oracletrees)
    samplesize = length(trainsamples)
    labelsize = targetsize(parser.labeler)

    i = 1
    testtrees = deepcopy(testsents)
    projectivize!(testtrees)
    testtrees = map(testtrees) do sent
        print("\r", i)
        i += 1
        s = State(sent, parser)
        while !isfinal(s)
            s = expandgold(s)[1]
        end
        s
    end

    info("OPTIMIZER: ", typeof(opt))
    info("BATCH SIZE: $batchsize LABEL SIZE: $labelsize")
    info("#BATCHES: $(div(samplesize, batchsize))")
    info("#SAMPLES: $(samplesize)")

    for i = 1:iter
        batch = sub(trainsamples, rand(1:samplesize, batchsize))
        preds = parser.labeler(map(b -> (b.w, b.t, b.l), batch))
        golds = map(s -> s.target, batch)
        correct = reduce(0, zip(argmax(preds.data, 1), golds)) do v, tup
            v + Int(tup[1] == tup[2])
        end
        accuracy = correct / batchsize
        loss = update!(opt, golds, preds)
        info("Epoch: $i\tLoss: $loss\tAcc: $accuracy")

        if i % evaliter == 0 && !isempty(testsents)
            info("Epoch: $i\tEvaluating labeler")
            @time res = assign!(parser.labeler, testtrees)
            res = convert(Vector{Vector{Int}}, res)
            uas, las = evaluate(parser, testtrees, res)

            info("Epoch: $i\tEvaluating pipeline")
            @time tree, labels = decode(LabelerFF, parser, testsents)
            uas, las = evaluate(parser, tree, labels)
            saver(parser, las) # decides whether to save based on LAS
        end
    end
end

# do pipeline job: unlabeled parser -> labeler
function decode{T}(::Type{LabelerFF}, parser::DepParser{T}, sents::Doc;
    batchsize=32, progbar=true)
    batches = 1:batchsize:length(sents)
    blabels = []
    res = []
    for k in batches
        batch = k:min(k+batchsize-1, endof(sents))
        ss = map(s -> State(s, parser), sub(sents, batch))
        map(setheads!, ss)
        parsegreedy!(parser, ss)
        labels = assign!(parser.labeler, ss)
        append!(res, ss)
        append!(blabels, labels)
    end
    res, blabels
end

max_stars_repo_name: hshindo/JukaiNLP

max_issues_repo_name: hshindo/JukaiNLP

max_forks_repo_name: hshindo/JukaiNLP


text:
module FresnelOptics

const VecI  = AbstractVector # Input  Vector
const VecO  = AbstractVector # Output Vector
const VecB  = AbstractVector # Buffer Vector
const VecIO = AbstractVector # In/Out Vector
const MatI  = AbstractMatrix # Input  Matrix
const MatO  = AbstractMatrix # Output Matrix
const MatB  = AbstractMatrix # Buffer Matrix
const MatIO = AbstractMatrix # In/Out Matrix

const 𝚷 = 2.0 * π

using LinearAlgebra.BLAS: gemv!

include("./utils/la.jl")
include("./fresnel.jl")
include("./lorentz.jl")
include("./dielectric.jl")
include("./layers.jl")

end # module