text:
"""
$(TYPEDEF)

A structure used to store information related to the bounds assigned to each
variable.

$(TYPEDFIELDS)
"""
Base.@kwdef struct VariableInfo{T<:AbstractFloat}
    "Is the variable integer valued?"
    is_integer::Bool                                            = false
    "Boolean indicating whether finite lower bound exists."
    has_lower_bound::Bool                                       = false
    "Boolean indicating whether finite upper bound exists."
    has_upper_bound::Bool                                       = false
    "Boolean indicating variable is fixed to a finite value."
    is_fixed::Bool                                              = false
    "Indicates that constraints have been set"
    has_constraints::Bool                                       = false
    "Lower bounds. May be -Inf."
    lower_bound::T                                              = typemin(T)
    "Upper bounds. May be Inf."
    upper_bound::T                                              = typemax(T)
end
is_integer(x::VariableInfo) = x.is_integer
has_lower_bound(x::VariableInfo) = x.has_lower_bound
has_upper_bound(x::VariableInfo) = x.has_upper_bound
lower_bound(x::VariableInfo{T}) where {T <: AbstractFloat} = x.lower_bound
upper_bound(x::VariableInfo{T}) where {T <: AbstractFloat} = x.upper_bound

is_fixed(x::VariableInfo) = x.is_fixed
function is_less_than(x::VariableInfo)
    flag = x.has_upper_bound
    flag &= !x.has_lower_bound
    return flag
end
function is_greater_than(x::VariableInfo)
    flag = x.has_lower_bound
    flag &= !x.has_upper_bound
    return flag
end
function is_interval(x::VariableInfo)
    flag = x.has_lower_bound
    flag &= x.has_upper_bound
    return flag
end
function is_int_interval(x::VariableInfo)
    flag = x.has_lower_bound
    flag &= x.has_upper_bound
    flag &= x.is_integer
    return flag
end
function is_real_interval(x::VariableInfo)
    flag = x.has_lower_bound
    flag &= x.has_upper_bound
    flag &= !x.is_integer
    return flag
end
function is_zero_one(x::VariableInfo{T}) where {T <: AbstractFloat}
    flag = iszero(x.lower_bound)
    flag &= isone(x.upper_bound)
    flag &= x.is_integer
    return flag
end

mid(x::VariableInfo{T}) where {T <: AbstractFloat} = 0.5*(upper_bound(x) + lower_bound(x))
diam(x::VariableInfo{T}) where {T <: AbstractFloat} = upper_bound(x) - lower_bound(x)
empty_variable_info(::Type{T}) where T = VariableInfo{T}(lower_bound = Inf,
                                                         upper_bound = -Inf)

Base.isempty(v::VariableInfo{T}) where {T <: AbstractFloat} = lower_bound(v) > upper_bound(v)
function check_isempty(l, u, b)
    flag = l <= u
    if b
        flag &= (l <= 0.0) | (u >= 1.0)
    end
    return !flag
end

function VariableInfo(::Type{T}, ::ZO) where {T <: AbstractFloat}
    return VariableInfo{T}(is_integer = true,
                           has_lower_bound = true,
                           has_upper_bound = true,
                           has_constraints = true,
                           lower_bound = zero(T),
                           upper_bound = one(T))
end

function VariableInfo(it::MOI.Interval{T}) where {T <: AbstractFloat}
    VariableInfo{T}(has_lower_bound = !isinf(it.lower),
                    has_upper_bound = !isinf(it.upper),
                    has_constraints = !isinf(it.lower) | !isinf(it.upper),
                    is_fixed = it.lower == it.upper,
                    lower_bound = it.lower,
                    upper_bound = it.upper)
end

function VariableInfo(v::VariableInfo{T}, ::ZO) where {T <: AbstractFloat}
    isempty(v) && (return v)
    l = max(zero(T), lower_bound(v))
    u = min(one(T), upper_bound(v))
    check_isempty(l, u, is_integer(v)) && return empty_variable_info(T)
    return VariableInfo{T}(is_integer = true,
                        has_lower_bound = true,
                        has_upper_bound = true,
                        has_constraints = true,
                        is_fixed = l == u,
                        lower_bound = l,
                        upper_bound = u)
end

function VariableInfo(v::VariableInfo{T}, it::MOI.Interval{T}) where {T <: AbstractFloat}
    isempty(v) && return v
    l = max(it.lower, lower_bound(v))
    u = min(it.upper, upper_bound(v))
    check_isempty(l, u, is_integer(v)) && return empty_variable_info(T)
    return VariableInfo(is_integer = is_integer(v),
                        has_lower_bound = !isinf(l),
                        has_upper_bound = !isinf(u),
                        has_constraints = !isinf(l) | !isinf(u),
                        is_fixed = l == u,
                        lower_bound = l,
                        upper_bound = u)
end

function VariableInfo(v::VariableInfo{T}, et::MOI.EqualTo{T}) where {T <: AbstractFloat}
    isempty(v) && return v
    l = max(et.value, lower_bound(v))
    u = min(et.value, upper_bound(v))
    check_isempty(l, u, is_integer(v)) && return empty_variable_info(T)
    out = VariableInfo(is_integer = is_integer(v),
    has_lower_bound = !isinf(l),
    has_upper_bound = !isinf(u),
    has_constraints = !isinf(l) | !isinf(u),
    is_fixed = true,
    lower_bound = l,
    upper_bound = u)
    return VariableInfo(is_integer = is_integer(v),
                        has_lower_bound = !isinf(l),
                        has_upper_bound = !isinf(u),
                        has_constraints = !isinf(l) | !isinf(u),
                        is_fixed = true,
                        lower_bound = l,
                        upper_bound = u)
end

function VariableInfo(v::VariableInfo{T}, gt::MOI.GreaterThan{T}) where {T <: AbstractFloat}
    isempty(v) && return v
    l = max(gt.lower, lower_bound(v))
    u = upper_bound(v)
    check_isempty(l, u, is_integer(v)) && return empty_variable_info(T)
    return VariableInfo(is_integer = is_integer(v),
                        has_lower_bound = !isinf(l),
                        has_upper_bound = !isinf(u),
                        has_constraints = !isinf(l) | !isinf(u),
                        is_fixed = l == u,
                        lower_bound = l,
                        upper_bound = u)
end
function VariableInfo(v::VariableInfo{T}, lt::MOI.LessThan{T}) where {T <: AbstractFloat}
    isempty(v) && return v
    l = lower_bound(v)
    u = min(lt.upper, upper_bound(v))
    check_isempty(l, u, is_integer(v)) && return empty_variable_info(T)
    return VariableInfo(is_integer = is_integer(v),
                        has_lower_bound = !isinf(l),
                        has_upper_bound = !isinf(u),
                        has_constraints = !isinf(l) | !isinf(u),
                        is_fixed = l == u,
                        lower_bound = l,
                        upper_bound = u)
end

function VariableInfo(v::VariableInfo{T}, s::MOI.Integer) where {T <: AbstractFloat}
    isempty(v) && return v
    l = lower_bound(v)
    u = upper_bound(v)
    check_isempty(l, u, true) && return empty_variable_info(T)
    return VariableInfo(is_integer = true,
                        has_lower_bound = has_lower_bound(v),
                        has_upper_bound = has_upper_bound(v),
                        has_constraints = v.has_constraints,
                        is_fixed = is_fixed(v),
                        lower_bound = l,
                        upper_bound = u)
end

ZO(v::VariableInfo)  = MOI.ZeroOne()
ET(v::VariableInfo{T}) where {T <: AbstractFloat}  = MOI.EqualTo{T}(v.lower_bound)
IT(v::VariableInfo{T}) where {T <: AbstractFloat}  = MOI.Interval{T}(v.lower_bound,v.upper_bound)
GT(v::VariableInfo{T}) where {T <: AbstractFloat}  = MOI.GreaterThan{T}(v.lower_bound)
LT(v::VariableInfo{T}) where {T <: AbstractFloat}  = MOI.LessThan{T}(v.upper_bound)
INT(v::VariableInfo{T}) where {T <: AbstractFloat} = MOI.Semiinteger{T}(v.lower_bound, v.upper_bound)

max_stars_repo_name: MatthewStuber/EAGO

max_issues_repo_name: MatthewStuber/EAGO

max_forks_repo_name: MatthewStuber/EAGO


text:
# test function from a debugging case in simScript.jl
# in julia/debug/bigSimulation/n6/julia
# seed 4545
# Claudia August 2015

if !isdefined(:individualtest) individualtest = false; end

if(individualtest)
    include("../src/types.jl")
    include("../src/functions.jl")
end

# to debug problem
tree="(3,(2,(((6,(5)#H9:0.91507):0.93066,(4,#H9:0.0):0.73688):0.0)#H7:1.79104::0.99498):0.11675,(1,#H7:0.04487::0.00502):0.4897);"
net0=readTopologyLevel1(tree);
#printEdges(net0)
net0.node[6].gammaz =1.0  #-5
net0.node[8].gammaz =0.067 #-7

q1 = Quartet(1,["6","5","4","1"],[0.5,0.4,0.1]);
qnet = extractQuartet!(net0,q1);
#printEdges(qnet)

identifyQuartet!(qnet)
qnet.which != 2 ? error("qnet which not correctly assigned") : nothing
#printEdges(qnet)

eliminateHybridization!(qnet)
qnet.which != 2 ? error("qnet which not correctly assigned") : nothing
qnet.numHybrids == 1 || error("qnet not correctly eliminated hyb")
qnet.hybrid[1].k == 4 || error("qnet not correclty identified")
qnet.hybrid[1].typeHyb == 5 || error("qnet not correclty identified")
qnet.hybrid[1].isBadDiamondI || error("qnet forgot it is a bad diamondI")

updateSplit!(qnet)
qnet.split != [-1,-1,-1,-1] ? error("qnet.split not correctly assigned") : nothing

updateFormula!(qnet)
qnet.formula != [-1,-1,-1] ? error("qnet.formula not correctly assigned") :nothing

calculateExpCF!(qnet)
qnet.expCF[2] < 0. || error("expCF not correctly calculated")

ext: jl


text:
module StockCorr

if VERSION >= v"0.7.0-DEV.3449"
    using LinearAlgebra
end

# The original, unoptimised code that simulates two correlated assets
function perf_stockcorr()
    ## Correlated asset information
    CurrentPrice = [78. 102.] # Initial Prices of the two stocks
    Corr = [1. 0.4; 0.4 1.]   # Correlation Matrix
    T = 500                   # Number of days to simulate = 2years = 500days
    n = 10000                 # Number of simulations
    dt = 1/250                # Time step (1year = 250days)
    Div = [0.01 0.01]         # Dividend
    Vol = [0.2 0.3]           # Volatility

    ## Market Information
    r = 0.03 # Risk-free rate

    ## Define storages
    SimulPriceA = zeros(T,n) # Simulated Price of Asset A
    SimulPriceA[1,:] .= CurrentPrice[1]
    SimulPriceB = zeros(T,n) # Simulated Price of Asset B
    SimulPriceB[1,:] .= CurrentPrice[2]

    ## Generating the paths of stock prices by Geometric Brownian Motion
    @static if VERSION <= v"0.7.0-DEV.5211"
        UpperTriangle = chol(Corr) # UpperTriangle Matrix by Cholesky decomposition
    else
        UpperTriangle = cholesky(Corr).U
    end

    for i = 1:n
       Wiener = randn(T-1,2)
       CorrWiener = Wiener*UpperTriangle
       for j = 2:T
          SimulPriceA[j,i] = SimulPriceA[j-1,i]*exp((r-Div[1]-Vol[1]^2/2)*dt+Vol[1]*sqrt(dt)*CorrWiener[j-1,1])
          SimulPriceB[j,i] = SimulPriceB[j-1,i]*exp((r-Div[2]-Vol[2]^2/2)*dt+Vol[2]*sqrt(dt)*CorrWiener[j-1,2])
       end
    end

    return (SimulPriceA, SimulPriceB)
end

end # module


text:
if VERSION >= v"0.5-"
    using Base.Test
else
    using BaseTestNext
    const Test = BaseTestNext
end

import Swifter: QueryResult, ResultInfo, App

@test 3 == QueryResult(ResultInfo(:any, 3), App(""),"",Dict())
@test 1.2 == QueryResult(ResultInfo(:any, 1.2), App(""),"",Dict())
@test true == QueryResult(ResultInfo(:any, true), App(""),"",Dict())
@test "Hello" == QueryResult(ResultInfo(:any, "Hello"), App(""),"",Dict())
@test AbstractString["A"] == QueryResult(ResultInfo(:any, AbstractString["A"]), App(""),"",Dict())
@test AbstractString["가"] == QueryResult(ResultInfo(:any, AbstractString["가"]), App(""),"",Dict())
@test Dict(1=>2) == QueryResult(ResultInfo(:any, Dict(1=>2)), App(""),"",Dict())

@test false == (3 === QueryResult(ResultInfo(:any, 3), App(""),"",Dict()))


text:
function myfuncB()
    println("Say hello from function B")
end

max_stars_repo_name: wigging/julia-computing

max_issues_repo_name: wigging/julia-computing

max_forks_repo_name: wigging/julia-computing


text:
module Assignment07

export normalizeDNA

# # uncomment the following line if you intend to use BioSequences types
# using BioSequences

"""
    normalizeDNA(::AbstractString)

Ensures that a sequence only contains valid bases
(or `'N'` for unknown bases).
Returns a String.
"""
function normalizeDNA(seq)
    seq = uppercase(string(seq))
    for base in seq
        # note: `N` indicates an unknown base
        occursin(base, "AGCTN") || error("invalid base $base")
    end
    return seq # change to `return LongDNASeq(seq)` if you want to try to use BioSequences types
end

# Your code here.
# Don't forget to export your functions!

end # module Assignment07


text:
# Interactions of the Julia and Python garbage collection:
#
#    * When we wrap a Python object in a Julia object jo (e.g. a PyDict),
#      we keep an explicit PyObject reference inside jo, whose finalizer
#      decrefs the Python object when it is called.
#
#    * When we wrap a Julia object jo inside a Python object po
#      (e.g a numpy array), we add jo to the pycall_gc dictionary,
#      keyed by a weak reference to po.  The Python weak reference
#      allows us to register a callback function that is called
#      when po is deallocated, and this callback function removes
#      jo from pycall_gc so that Julia can garbage-collect it.

pycall_gc = Dict{PyPtr,Any}()

function weakref_callback(callback::PyPtr, wo::PyPtr)
    global pycall_gc
    try
        delete!(pycall_gc::Dict{PyPtr,Any}, wo)
        # not sure what to do if there is an exception here
    finally
        ccall((@pysym :Py_DecRef), Void, (PyPtr,), wo)
    end
    ccall((@pysym :Py_IncRef), Void, (PyPtr,), pynothing::PyPtr)
    return pynothing::PyPtr
end

weakref_callback_obj = PyObject() # weakref_callback Python method

function pygc_finalize()
    global pycall_gc
    global weakref_callback_obj
    pydecref(weakref_callback_obj)
    pycall_gc::Dict{PyPtr,Any} = Dict{PyPtr,Any}()
end

# "embed" a reference to jo in po, using the weak-reference mechanism
function pyembed(po::PyObject, jo::Any)
    global pycall_gc
    global weakref_callback_obj
    if (weakref_callback_obj::PyObject).o == C_NULL
        weakref_callback_obj::PyObject = pymethod(weakref_callback,
                                                  "weakref_callback",
                                                  METH_O)
    end
    wo = @pycheckn ccall((@pysym :PyWeakref_NewRef), PyPtr, (PyPtr,PyPtr),
                         po, weakref_callback_obj)
    (pycall_gc::Dict{PyPtr,Any})[wo] = jo
    return po
end


text:
# Modified scheme of Bassi and Rebay (BR2)

"""
  Penalty for BR2
"""
struct BR2Penalty{Tsol, Tres} <: AbstractDiffusionPenalty

  # primal method
  delta_w_n::Matrix{Tsol}
  qL::Array{Tres, 3}
  qR::Array{Tres, 3}
  t1L::Array{Tres, 3}
  t1R::Array{Tres, 3}
  t2L::Array{Tres, 3}
  t2R::Array{Tres, 3}

  # diff method
  t1_dotL::Array{Tres, 5}
  t1_dotR::Array{Tres, 5}

  # the first L/R indicates if this is element kappa or nu
  # the second L/R indices if the derivative is wrt qL or qR
  t2LL::Array{Tres, 5}
  t2LR::Array{Tres, 5}
  t2RL::Array{Tres, 5}
  t2RR::Array{Tres, 5}

  t3LL::Array{Tres, 5}
  t3LR::Array{Tres, 5}
  t3RL::Array{Tres, 5}
  t3RR::Array{Tres, 5}

  # re-use t1_dotL,R
  t4RL::Array{Tres, 5}
  t4RR::Array{Tres, 5}

  function BR2Penalty{Tsol, Tres}(mesh, sbp, opts) where {Tsol, Tres}

    @unpack mesh numDofPerNode numNodesPerFace numNodesPerElement dim

    delta_w_n = zeros(Tsol, numDofPerNode, numNodesPerFace)
    qL = zeros(Tres, numDofPerNode, numNodesPerElement, dim)
    qR = zeros(Tres, numDofPerNode, numNodesPerElement, dim)
    t1L = zeros(Tres, numDofPerNode, numNodesPerElement, dim)
    t1R = zeros(Tres, numDofPerNode, numNodesPerElement, dim)
    t2L = zeros(Tres, numDofPerNode, numNodesPerFace, dim)
    t2R = zeros(Tres, numDofPerNode, numNodesPerFace, dim)

    t1_dotL = zeros(Tres, numDofPerNode, numDofPerNode, dim, numNodesPerFace,
                   numNodesPerElement)
    t1_dotR = zeros(Tres, numDofPerNode, numDofPerNode, dim, numNodesPerFace,
                     numNodesPerElement)

    t2LL = zeros(Tres, numDofPerNode, numDofPerNode, dim, numNodesPerElement,
                      numNodesPerElement)
    t2LR = zeros(Tres, numDofPerNode, numDofPerNode, dim, numNodesPerElement,
                      numNodesPerElement)

    t2RL = zeros(Tres, numDofPerNode, numDofPerNode, dim, numNodesPerElement,
                      numNodesPerElement)
    t2RR = zeros(Tres, numDofPerNode, numDofPerNode, dim, numNodesPerElement,
                      numNodesPerElement)

    t3LL = zeros(Tres, numDofPerNode, numDofPerNode, dim, numNodesPerElement,
                       numNodesPerElement)
    t3LR = zeros(Tres, numDofPerNode, numDofPerNode, dim, numNodesPerElement,
                       numNodesPerElement)

    t3RL = zeros(Tres, numDofPerNode, numDofPerNode, dim, numNodesPerElement,
                       numNodesPerElement)
    t3RR = zeros(Tres, numDofPerNode, numDofPerNode, dim, numNodesPerElement,
                       numNodesPerElement)

    t4RL = zeros(Tres, numDofPerNode, numDofPerNode, dim, numNodesPerFace,
                     numNodesPerElement)
    t4RR = zeros(Tres, numDofPerNode, numDofPerNode, dim, numNodesPerFace,
                       numNodesPerElement)

    return new(delta_w_n, qL, qR, t1L, t1R, t2L, t2R,
               t1_dotL, t1_dotR, t2LL, t2LR, t2RL, t2RR, t3LL, t3LR, t3RL, t3RR,
               t4RL, t4RR)
  end

end

"""
  Applies the penalty for the SBP-SAT generalization of the modified scheme
  of Bassi and Rebay (BR2).
"""
function applyPenalty(penalty::BR2Penalty{Tsol, Tres}, sbp, params::AbstractParamType,
                      sbpface,
                      diffusion::AbstractDiffusion, iface::Interface,
                      delta_w::AbstractMatrix{Tsol}, theta::AbstractMatrix{Tres},
                      qL_el::AbstractMatrix, qR_el::AbstractMatrix,
                      wL::AbstractMatrix, wR::AbstractMatrix,
                      coordsL::AbstractMatrix, coordsR::AbstractMatrix,
                      nrm_face::AbstractMatrix,
                      alphas::AbstractVector,
                      dxidxL::Abstract3DArray, dxidxR::Abstract3DArray,
                      jacL::AbstractVector, jacR::AbstractVector,
                      res1L::AbstractMatrix, res1R::AbstractMatrix,
                      res2L::AbstractMatrix, res2R::AbstractMatrix
                     ) where {Tsol, Tres}

  fill!(res1L, 0); fill!(res1R, 0)
  fill!(res2L, 0); fill!(res2R, 0)
  numDofPerNode, numNodesPerFace = size(delta_w)
  numNodesPerElement = length(jacL)
  dim = size(nrm_face, 1)

  @unpack penalty delta_w_n qL qR t1L t1R t2L t2R
  fill!(qL, 0); fill!(qR, 0)

  #--------------------------
  # apply T1
  # multiply by normal vector, then R^T B
  #alpha_g = 1/(dim + 1)  # = 1/number of faces of a simplex
  @simd for d1=1:dim
    @simd for j=1:numNodesPerFace
      @simd for k=1:numDofPerNode
        delta_w_n[k, j] = 0.25*delta_w[k, j]*nrm_face[d1, j]
      end
    end

    qL_d = sview(qL, :, :, d1); qR_d = sview(qR, :, :, d1)
    interiorFaceIntegrate!(sbpface, iface, delta_w_n, qL_d, qR_d)

    # interiorFaceIntegrates -= the second output
    @simd for j=1:numNodesPerElement
      @simd for k=1:numDofPerNode
        qR_d[k, j] = -qR_d[k, j]
      end
    end
  end

  # apply Lambda matrix
  applyDiffusionTensor(diffusion, sbp, params, qL_el, wL, coordsL, dxidxL, jacL, iface.elementL, sbpface, iface.faceL,
                       qL, t1L)
  applyDiffusionTensor(diffusion, sbp, params, qR_el, wR, coordsR, dxidxR, jacR,
                       iface.elementR, sbpface, iface.faceR,
                       qR, t1R)

  # apply inverse mass matrix, then apply B*Nx*R*t2L_x + B*Ny*R*t2L_y
  @simd for d1=1:dim
    @simd for j=1:numNodesPerElement
      facL = (1/alphas[1])*jacL[j]/sbp.w[j]
      facR = (1/alphas[2])*jacR[j]/sbp.w[j]
      #facL = (3)*jacL[j]/sbp.w[j]
      #facR = (3)*jacR[j]/sbp.w[j]

      @simd for k=1:numDofPerNode
        t1L[k, j, d1] *= facL
        t1R[k, j, d1] *= facR
      end
    end

    #TODO: make t2L 2D rather than 3D
    t1L_d = ro_sview(t1L, :, :, d1); t1R_d = ro_sview(t1R, :, :, d1)
    t2L_d = sview(t2L, :, :, d1);    t2R_d = sview(t2R, :, :, d1)
    interiorFaceInterpolate!(sbpface, iface, t1L_d, t1R_d, t2L_d, t2R_d)

    @simd for j=1:numNodesPerFace
      @simd for k=1:numDofPerNode
        val = sbpface.wface[j]*nrm_face[d1, j]*(t2L_d[k, j] + t2R_d[k, j])
        res1L[k, j] += val
        res1R[k, j] -= val  # delta_w is reversed for elementR
      end
    end

  end  # end d1

  #----------------------------
  # apply T2 and T3
  @simd for j=1:numNodesPerFace
    @simd for k=1:numDofPerNode
      val1 =  0.5*sbpface.wface[j]*theta[k, j]
      val2 = -0.5*sbpface.wface[j]*delta_w[k, j]
      res1L[k, j] += val1
      res1R[k, j] += val1
      res2L[k, j] += val2
      res2R[k, j] -= val2  # delta_w is reversed for elementR
    end
  end

  return nothing
end

"""
  Applies the Dirichlet penalty matrix T_gammaD.

  **Inputs**

   * penalty: [`AbstractDiffusionPenalty`](@ref)
   * sbp
   * sbpface
   * diffusion: [`AbstractDiffusion`](@ref)
   * bndry: `Boundary` object
   * delta_w: values to multiply the penalty against, `numDofPerNode` x
              `numNodesPerFace`
   * wL: entropy variables for the element, `numDofPerNode` x `numNodesPerElement`
   * nrm_face: (scaled) normal vector at each face node, `dim` x
               `numNodesPerFace`
   * alpha: the alpha_gk parameter (Number)
   * jacL: mapping jacobian determinant, `numNodesPerElement`

  **Inputs/Outputs**

   * res1L: `numDofPerNode` x `numNodesPerFace` array to overwrite with
            the result

"""
#=
function applyDirichletPenalty(penalty::BR2Penalty{Tsol, Tres}, sbp, sbpface,
                      diffusion::AbstractDiffusion, bndry::Boundary,
                      delta_w::AbstractMatrix{Tsol},
                      wL::AbstractMatrix,
                      nrm_face::AbstractMatrix,
                      alpha::Number,
                      jacL::AbstractVector,
                      res1L::AbstractMatrix,
                     ) where {Tsol, Tres}
=#
function applyDirichletPenalty(penalty::BR2Penalty{Tsol, Tres}, sbp,
                      params::AbstractParamType, sbpface,
                      diffusion::AbstractDiffusion, bndry::Boundary,
                      delta_w::AbstractMatrix,
                      uL::AbstractMatrix,
                      wL::AbstractMatrix,
                      coordsL::AbstractMatrix,
                      nrm_face::AbstractMatrix,
                      alpha::Number,
                      dxidxL::Abstract3DArray,
                      jacL::AbstractVector,
                      res1L::AbstractMatrix)  where {Tsol, Tres}

  fill!(res1L, 0)
  numDofPerNode, numNodesPerFace = size(delta_w)
  numNodesPerElement = length(jacL)
  dim = size(nrm_face, 1)

  @unpack penalty delta_w_n qL qR t1L t1R t2L t2R
  fill!(qL, 0); fill!(qR, 0)

  # multiply by normal vector, then R^T B
  #alpha_g = 1/(dim + 1)  # = 1/number of faces of a simplex
  @simd for d1=1:dim
    @simd for j=1:numNodesPerFace
      @simd for k=1:numDofPerNode
        delta_w_n[k, j] = delta_w[k, j]*nrm_face[d1, j]
      end
    end

    qL_d = sview(qL, :, :, d1); qR_d = sview(qR, :, :, d1)
    boundaryFaceIntegrate!(sbpface, bndry.face, delta_w_n, qL_d)
  end

  # apply Lambda matrix
  applyDiffusionTensor(diffusion, sbp, params, uL, wL, coordsL, dxidxL, jacL,
                       bndry.element, sbpface, bndry.face,
                       qL, t1L)

  # apply inverse mass matrix, then apply B*Nx*R*t2L_x + B*Ny*R*t2L_y
  @simd for d1=1:dim
    @simd for j=1:numNodesPerElement
      facL = (1/alpha)*jacL[j]/sbp.w[j]
      #facL = (3)*jacL[j]/sbp.w[j]
      #facR = (3)*jacR[j]/sbp.w[j]

      @simd for k=1:numDofPerNode
        t1L[k, j, d1] *= facL
      end
    end

    #TODO: make t2L 2D rather than 3D
    t1L_d = ro_sview(t1L, :, :, d1)
    t2L_d = sview(t2L, :, :, d1)
    boundaryFaceInterpolate!(sbpface, bndry.face, t1L_d, t2L_d)

    @simd for j=1:numNodesPerFace
      @simd for k=1:numDofPerNode
        val = sbpface.wface[j]*nrm_face[d1, j]*t2L_d[k, j]
        res1L[k, j] += val
      end
    end

  end  # end d1

  return nothing
end

"""
  Differentiated version of BR2 penalty
"""
function applyPenalty_diff(penalty::BR2Penalty{Tsol, Tres}, sbp,
                      params::AbstractParamType, sbpface,
                      diffusion::AbstractDiffusion, iface::Interface,
                      delta_w::AbstractMatrix{Tsol},
                      delta_w_dotL::Abstract4DArray, delta_w_dotR::Abstract4DArray,
                      theta::AbstractMatrix{Tres},
                      theta_dotL::Abstract4DArray, theta_dotR::Abstract4DArray,
                      qL_el::AbstractMatrix, qR_el::AbstractMatrix,
                      wL::AbstractMatrix, wR::AbstractMatrix,
                      coordsL::AbstractMatrix, coordsR::AbstractMatrix,
                      nrm_face::AbstractMatrix,
                      alphas::AbstractVector,
                      dxidxL::Abstract3DArray, dxidxR::Abstract3DArray,
                      jacL::AbstractVector, jacR::AbstractVector,
                      res1L_dotL::Abstract4DArray, res1L_dotR::Abstract4DArray,
                      res1R_dotL::Abstract4DArray, res1R_dotR::Abstract4DArray,
                      res2L::AbstractMatrix, res2R::AbstractMatrix,
                      res2L_dotL::Abstract4DArray, res2L_dotR::Abstract4DArray,
                      res2R_dotL::Abstract4DArray, res2R_dotR::Abstract4DArray,
                     ) where {Tsol, Tres}

  numDofPerNode, numNodesPerFace = size(delta_w)
  numNodesPerElement = length(jacL)
  dim = size(nrm_face, 1)
  add = SummationByParts.Add()

  # nodes in the stencil of R
  nodesL = sview(sbpface.perm, :, iface.faceL)
  nodesR = sview(sbpface.perm, :, iface.faceR)

  # the first L/R indicates if this is element kappa or nu
  # the second L/R indices if the derivative is wrt qL or qR
  @unpack penalty delta_w_n qL qR t1_dotL t1_dotR t2LL t2LR t2RL t2RR
  @unpack penalty t3LL t3LR t3RL t3RR t4RL t4RR
  t4LL = t1_dotL; t4LR = t1_dotR

  fill!(qL, 0); fill!(qR, 0)

  #--------------------------
  # T1

  # compute primal quantity needed later
  # multiply by normal vector, then R^T B
  #alpha_g = 1/(dim + 1)  # = 1/number of faces of a simplex
  @simd for d1=1:dim
    @simd for j=1:numNodesPerFace
      @simd for k=1:numDofPerNode
        delta_w_n[k, j] = 0.25*delta_w[k, j]*nrm_face[d1, j]
      end
    end

    qL_d = sview(qL, :, :, d1); qR_d = sview(qR, :, :, d1)
    interiorFaceIntegrate!(sbpface, iface, delta_w_n, qL_d, qR_d)

    # interiorFaceIntegrates -= the second output
    @simd for j=1:numNodesPerElement
      @simd for k=1:numDofPerNode
        qR_d[k, j] = -qR_d[k, j]
      end
    end
  end

  # compute derivative quantity
  # apply B and [Nx, Ny]
  @simd for q=1:numNodesPerElement
    @simd for p=1:numNodesPerFace
      @simd for d=1:dim
        fac = 0.25*sbpface.wface[p]*nrm_face[d, p]
        @simd for j=1:numDofPerNode
          @simd for i=1:numDofPerNode
            t1_dotL[i, j, d, p, q] = fac*delta_w_dotL[i, j, p, q]
            t1_dotR[i, j, d, p, q] = fac*delta_w_dotR[i, j, p, q]
          end
        end
      end
    end
  end

  # interpolate back to volume nodes
  fill!(t2LL, 0); fill!(t2RL, 0); fill!(t2LR, 0); fill!(t2RR, 0)
  interiorFaceIntegrate_jac!(sbpface, iface, t1_dotL, t1_dotL, t2LL, t2RL, add,
                             add, false)
  interiorFaceIntegrate_jac!(sbpface, iface, t1_dotR, t1_dotR, t2LR, t2RR, add,
                             add, false)

  # apply diffusion tensor
  applyDiffusionTensor_diff(diffusion, sbp, params, qL_el, wL, coordsL, dxidxL,
                            jacL, iface.elementL, sbpface,
                            iface.faceL, qL, t2LL, t2LR, t3LL, t3LR)
  # reverse RR and RL here because Lambda is a function of qR, not qL
  applyDiffusionTensor_diff(diffusion, sbp, params, qR_el, wR, coordsR, dxidxR,
                            jacR, iface.elementR, sbpface,
                            iface.faceR, qR, t2RR, t2RL, t3RR, t3RL)

  # apply inverse mass matrix (only nodes needed for interpolation later)
  @simd for q=1:numNodesPerElement
    # elementL
    @simd for p in nodesL
      facL = (1/alphas[1])*jacL[p]/sbp.w[p]
      @simd for d=1:dim
        @simd for j=1:numDofPerNode
          @simd for i=1:numDofPerNode
            t3LL[i, j, d, p, q] *= facL
            t3LR[i, j, d, p, q] *= facL
          end
        end
      end
    end
    # elementR
    @simd for p in nodesR
      facR = (1/alphas[2])*jacR[p]/sbp.w[p]
      @simd for d=1:dim
        @simd for j=1:numDofPerNode
          @simd for i=1:numDofPerNode
            t3RL[i, j, d, p, q] *= facR
            t3RR[i, j, d, p, q] *= facR
          end
        end
      end
    end
  end  # end q

  # interpolate back to the face
  fill!(t4LL, 0); fill!(t4RL, 0); fill!(t4LR, 0); fill!(t4RR, 0)
  interiorFaceInterpolate_jac!(sbpface, iface, t3LL, t3RL, t4LL, t4RL)
  interiorFaceInterpolate_jac!(sbpface, iface, t3LR, t3RR, t4LR, t4RR)

  # apply B * [Nx, Ny] and sum
#  fill!(res1L_dotL, 0); fill!(res1L_dotR, 0)  #TODO: tile this
#  fill!(res1R_dotL, 0); fill!(res1R_dotR, 0)
  @simd for q=1:numNodesPerElement
    @simd for p=1:numNodesPerFace
      @simd for j=1:numDofPerNode
        @simd for i=1:numDofPerNode
          res1L_dotL[i, j, p, q] = 0
          res1L_dotR[i, j, p, q] = 0
          res1R_dotL[i, j, p, q] = 0
          res1R_dotR[i, j, p, q] = 0
        end
      end

      @simd for d=1:dim
        fac = sbpface.wface[p]*nrm_face[d, p]
        @simd for j=1:numDofPerNode
          @simd for i=1:numDofPerNode
            valL = fac*(t4LL[i, j, d, p, q] + t4RL[i, j, d, p, q])
            valR = fac*(t4LR[i, j, d, p, q] + t4RR[i, j, d, p, q])
            #println("t2L_d_dot = ", real(t4LL[i, j, d, p, q]))
            #println("t2R_d_dot = ", real(t4RL[i, j, d, p, q]))
            #println("val = ", real(valL))
            res1L_dotL[i, j, p, q] +=  valL
            res1L_dotR[i, j, p, q] +=  valR
            res1R_dotL[i, j, p, q] += -valL  # delta_w is reversed for elementR
            res1R_dotR[i, j, p, q] += -valR
          end
        end
      end
    end
  end

  # T2 and T3
  @simd for q=1:numNodesPerElement
    @simd for p=1:numNodesPerFace
      @simd for j=1:numDofPerNode
        @simd for i=1:numDofPerNode
          res2L_dotL[i, j, p, q] = 0
          res2L_dotR[i, j, p, q] = 0
          res2R_dotL[i, j, p, q] = 0
          res2R_dotR[i, j, p, q] = 0
        end
      end
      fac1 =  0.5*sbpface.wface[p]
      fac2 = -0.5*sbpface.wface[p]
      @simd for j=1:numDofPerNode
        @simd for i=1:numDofPerNode
          res1L_dotL[i, j, p, q] += fac1*theta_dotL[i, j, p, q]
          res1L_dotR[i, j, p, q] += fac1*theta_dotR[i, j, p, q]
          res1R_dotL[i, j, p, q] += fac1*theta_dotL[i, j, p, q]
          res1R_dotR[i, j, p, q] += fac1*theta_dotR[i, j, p, q]

          res2L_dotL[i, j, p, q] += fac2*delta_w_dotL[i, j, p, q]
          res2L_dotR[i, j, p, q] += fac2*delta_w_dotR[i, j, p, q]
          res2R_dotL[i, j, p, q] -= fac2*delta_w_dotL[i, j, p, q]
          res2R_dotR[i, j, p, q] -= fac2*delta_w_dotR[i, j, p, q]
        end
      end
    end
  end

  # also need res2 as output
  fill!(res2L, 0); fill!(res2R, 0)
  @simd for j=1:numNodesPerFace
    @simd for k=1:numDofPerNode
      val2 = -0.5*sbpface.wface[j]*delta_w[k, j]
      res2L[k, j] =  val2
      res2R[k, j] = -val2  # delta_w is reversed for elementR
    end
  end

  return false
end

"""
  Differentiated version of BR2 penalty
"""
function applyDirichletPenalty_diff(penalty::BR2Penalty{Tsol, Tres}, sbp,
                      params::AbstractParamType,
                      sbpface, diffusion::AbstractDiffusion, bndry::Boundary,
                      delta_w::AbstractMatrix{Tsol},
                      delta_w_dot::Abstract4DArray,
                      qL_el::AbstractMatrix,
                      wL::AbstractMatrix,
                      coordsL::AbstractMatrix,
                      nrm_face::AbstractMatrix,
                      alpha::Number,
                      dxidxL::Abstract3DArray,
                      jacL::AbstractVector,
                      res_dot::Abstract4DArray,
                     ) where {Tsol, Tres}

  numDofPerNode, numNodesPerFace = size(delta_w)
  numNodesPerElement = length(jacL)
  dim = size(nrm_face, 1)
  add = SummationByParts.Add()

  # nodes in the stencil of R
  nodesL = sview(sbpface.perm, :, bndry.face)

  # the first L/R indicates if this is element kappa or nu
  # the second L/R indices if the derivative is wrt qL or qR
  @unpack penalty delta_w_n qL t1_dotL t2LL t3LL
  t4LL = t1_dotL

  fill!(qL, 0);

  # compute primal quantity needed later
  # multiply by normal vector, then R^T B
  #alpha_g = 1/(dim + 1)  # = 1/number of faces of a simplex
  @simd for d1=1:dim
    @simd for j=1:numNodesPerFace
      @simd for k=1:numDofPerNode
        delta_w_n[k, j] = delta_w[k, j]*nrm_face[d1, j]
      end
    end

    qL_d = sview(qL, :, :, d1)
    boundaryFaceIntegrate!(sbpface, bndry.face, delta_w_n, qL_d)
  end

  # compute derivative quantity
  # apply B and [Nx, Ny]
  @simd for q=1:numNodesPerElement
    @simd for p=1:numNodesPerFace
      @simd for d=1:dim
        fac = sbpface.wface[p]*nrm_face[d, p]
        @simd for j=1:numDofPerNode
          @simd for i=1:numDofPerNode
            t1_dotL[i, j, d, p, q] = fac*delta_w_dot[i, j, p, q]
          end
        end
      end
    end
  end

  # interpolate back to volume nodes
  fill!(t2LL, 0)
  boundaryFaceIntegrate_jac!(sbpface, bndry.face, t1_dotL, t2LL, include_quadrature=false)

  # apply diffusion tensor
  applyDiffusionTensor_diff(diffusion, sbp, params, qL_el, wL, coordsL, dxidxL,
                            jacL, bndry.element, sbpface,
                            bndry.face, qL, t2LL, t3LL)

  # apply inverse mass matrix (only nodes needed for interpolation later)
  @simd for q=1:numNodesPerElement
    # elementL
    @simd for p in nodesL
      facL = (1/alpha)*jacL[p]/sbp.w[p]
      @simd for d=1:dim
        @simd for j=1:numDofPerNode
          @simd for i=1:numDofPerNode
            t3LL[i, j, d, p, q] *= facL
          end
        end
      end
    end
  end  # end q

  # interpolate back to the face
  fill!(t4LL, 0)
  boundaryFaceInterpolate_jac!(sbpface, bndry.face, t3LL, t4LL)

  # apply B * [Nx, Ny] and sum
  @simd for q=1:numNodesPerElement
    @simd for p=1:numNodesPerFace
      @simd for j=1:numDofPerNode
        @simd for i=1:numDofPerNode
          res_dot[i, j, p, q] = 0
        end
      end

      @simd for d=1:dim
        fac = sbpface.wface[p]*nrm_face[d, p]
        @simd for j=1:numDofPerNode
          @simd for i=1:numDofPerNode
            res_dot[i, j, p, q] +=  fac*t4LL[i, j, d, p, q]
          end
        end
      end
    end
  end

  return nothing
end


text:
module SpreadBuilder

__precompile__(true)

using ..GreekAnalytical
using ..PositionTypes: Position
export portfolio_greeks, portfolio_PnL

function portfolio_greeks(portf::AbstractVector{Position}, underlying_price::T, r::T, b::T) where {T<:Real}
    ref_maturity = minimum(x->x.DTE, portf) / 365.0

    portf_delta = 0.0
    portf_gamma = 0.0
    portf_vega = 0.0
    portf_theta = 0.0
    portf_vanna = 0.0
    portf_volga = 0.0
    portf_charm = 0.0
    portf_color = 0.0

    for pos in portf
        if pos.holding isa Options
            option = pos.holding
            quantity = pos.quantity
            portf_delta += blackScholesDelta(option.isCall, underlying_price, option.strike, to_τ(option), r, b, get_ivol(option)) * quantity
            portf_gamma += blackScholesGamma(option.isCall, underlying_price, option.strike, to_τ(option), r, b, get_ivol(option)) * quantity
            portf_theta += blackScholesTheta(option.isCall, underlying_price, option.strike, to_τ(option), r, b, get_ivol(option)) * quantity
            portf_vega += sqrt(to_τ(option) / ref_maturity) * blackScholesVega(option.isCall, underlying_price, option.strike, to_τ(option), r, b, get_ivol(option)) * quantity
            portf_vanna += sqrt(to_τ(option) / ref_maturity) * blackScholesVanna(option.isCall, underlying_price, option.strike, to_τ(option), r, b, get_ivol(option)) * quantity
            portf_volga += sqrt(to_τ(option) / ref_maturity) * blackScholesVolga(option.isCall, underlying_price, option.strike, to_τ(option), r, b, get_ivol(option)) * quantity
            portf_charm += blackScholesCharm(option.isCall, underlying_price, option.strike, to_τ(option), r, b, get_ivol(option)) * quantity
            portf_color += blackScholesColor(option.isCall, underlying_price, option.strike, to_τ(option), r, b, get_ivol(option)) * quantity
        else
            quantity = pos.quantity
            portf_delta += pos.quantity
        end
    end
end

function portfolio_PnL(portf::AbstractVector{Position}, initial_price::T, underlying_price::T, r::T, b::T) where {T<:Real}
    pnl = 0.0

    for pos in portf
        if pos.holding isa Options
            option = pos.holding
            quantity = pos.quantity
            pnl += (blackScholesValue(option.isCall, underlying_price, option.strike, to_τ(option), r, b, get_ivol(option)) - blackScholesValue(option.isCall, initial_price, option.strike, to_τ(option), r, b, get_ivol(option))) * quantity
        else
            pnl += underlying_price - initial_price
        end
    end
end

end

max_stars_repo_name: Stochastic-Adventure/BlackScholesGreeks

max_issues_repo_name: Stochastic-Adventure/BlackScholesGreeks

max_forks_repo_name: Stochastic-Adventure/BlackScholesGreeks


text:
module ProjectEulerUtil

# Write your package code here.
include("proper_divisors.jl")
include("abundant-numbers.jl")
include("inverse_digits.jl")
include("get_problem_title.jl")

end


text:
export project_subspace!

"""
project a vector x onto the subspace spanned by the columns of A.
This is a simple implementation that solves linear systems on the fly. We can
add options to factor A' * A once and reuse if memory permits. Iterative factorization
free methods are also an option. Orthogonalize subspace beforehand to obtain an A
such that A' * A = I to avoid linear system solves all together in this stage.
"""
function project_subspace!(
                          x     ::Vector{TF},
                          A     ::Union{Array{Any},Array{TF,2}},
                          orth  ::Bool
                          ) where {TF<:Real}
  if orth == true
     x .= A*(A'*x)
  else
    x .= A*((A'*A)\(A'*x))
  end

end

function project_subspace!(
                          x     ::Array{TF,2},
                          A     ::Union{Array{Any},Array{TF,2}},
                          orth  ::Bool,
                          mode  ::Tuple{String,String}
                          ) where {TF<:Real}

(n1,n2) = deepcopy(size(x))

if mode[1] == "slice"
  error("mode[1] for project_subspace! must be: fiber")
end

if mode[2] == "x" #project each x[:,i] onto the subspace, so each column of x
  if orth == true
     x .= A*(A'*x)
  else
    x .= A*((A'*A)\(A'*x))
  end
elseif mode[2] == "z" #project each x[i,:] onto the subspace, so each row of x
  if orth == true
     x .= (A*(A'*x'))'
  else
    x .= (A*((A'*A)\(A'*x')))'
  end
else
  error("mode[2] for project_subspace! with 2D array input must be: x, or z")
end

x = vec(x)
end

"""
we want to project each slice of a tensor onto the subspace spanned by the columns of A
x is a 3D array and we need to permute and reshape it such that we project the
correct slices. Every slice needs to become a column of the matrix (x) we projection onto the subspace in column-wise sense
"""
function project_subspace!(
                          x     ::Array{TF,3},
                          A     ::Union{Array{Any},Array{TF,2}},
                          orth  ::Bool,
                          mode  ::Tuple{String,String}
                          ) where {TF<:Real}

  (n1,n2,n3) = deepcopy(size(x))

if mode[1] == "slice"
  #permute and reshape
  if mode[2] == "x"
    x = permutedims(x,[2;3;1]);
    x = reshape(x,n2*n3,n1)
  elseif mode[2] == "y"
    x = permutedims(x,[1;3;2]);
    x = reshape(x,n1*n3,n2)
  elseif mode[2] == "z"
    x = reshape(x,n1*n2,n3) #there are n3 slices of dimension n1*n2
  end

  #project every column of matrix x (a slice of the original tensor)
  if orth == true
     x .= A*(A'*x)
  else
    x .= A*((A'*A)\(A'*x))
  end

  #permute back and vectorize
  if mode[2] == "x"
    x = reshape(x,n2,n3,n1)
    x = permutedims(x,[3;1;2]);
  elseif mode[2] == "y"
    x = reshape(x,n1,n3,n2)
    x = permutedims(x,[1;3;2]);
  elseif mode[2] == "z"
    x = reshape(x,n1,n2,n3) #there are n3 slices of dimension n1*n2
  end

else
  error("for 3D models, the mode of application for project_subspace! needs to be (slice,x) or (slice,y) or (slice,z)")

end
x = vec(x)
end

max_stars_repo_name: slimgroup/SetIntersectionProjection

max_issues_repo_name: slimgroup/SetIntersectionProjection

max_forks_repo_name: slimgroup/SetIntersectionProjection


text:
module TiramisuPretag

using ..ParallelAccelerator

import CompilerTools

export process_node

type State
  forCount::Integer
  ifCount::Integer
  currentFor::Integer
end

function process_node(node::Expr, state::Union{State, Void}, top_level_number, is_top_level, read)
  if state === nothing
    state = State(0, 0, -1)
  end

  parentFor = nothing

  if node.head === :for
    parentFor = state.currentFor
    state.currentFor = state.forCount
    state.forCount += 1
  elseif node.head === :if
    state.ifCount += 1
  end

  children = Vector{Any}()
  for arg in node.args
    if isa(arg, Expr)
      if arg.head === :for
        push!(children, Expr(:meta, :forloop, state.forCount))
      elseif arg.head === :continue
        push!(children, Expr(:meta, :continue, state.currentFor))
      elseif arg.head === :if
        push!(children, Expr(:meta, :if, state.ifCount))
        if 3 <= length(arg.args)  # If there is an else
          index = 1
          elseblock = arg.args[3].args
          if isa(elseblock[1], Expr) && elseblock[1].head === :line
            index = 2 # skip line info
          end
          insert!(elseblock, index, Expr(:meta, :else, state.ifCount))
          push!(elseblock, Expr(:meta, :endif, state.ifCount))
        else
          push!(arg.args[2].args, Expr(:meta, :endif, state.ifCount))
        end
      end
    end
    push!(children, CompilerTools.AstWalker.AstWalk(arg, process_node, state))
  end

  node.args = children

  if node.head === :for
    state.currentFor = parentFor
  end

  return node
end

function process_node(node::Any, state::Union{State, Void}, top_level_number, is_top_level, read)
  CompilerTools.AstWalker.ASTWALK_RECURSE
end

end

max_stars_repo_name: malekbr/julia-to-tiramisu

max_issues_repo_name: malekbr/julia-to-tiramisu

max_forks_repo_name: malekbr/julia-to-tiramisu


text:
using GenericInstruments
GI = GenericInstruments
####################################
# Instantiate obj
niFgen1 = GI.AWG.INSTR(:NI5421,"PXI1Slot2")
# Init
GI.AWG.connect!(niFgen1)
x=GI.AWG.get_idn(niFgen1)
GI.AWG.disconnect!(niFgen1)
# pending...


text:
@testset "smarts.smiles" begin

@testset "chain" begin
    nullmol = parse(SMILES, "")
    @test nodecount(nullmol) == 0
    @test edgecount(nullmol) == 0

    methane = parse(SMILES, "C")
    @test nodeattr(methane, 1).symbol == :C
    @test edgecount(methane) == 0

    EtOH = parse(SMILES, "CCO")
    @test nodeattr(EtOH, 3).symbol == :O

    EtOH2 = parse(SMILES, "[CH3][CH2][OH]")
    @test nodecount(EtOH2) == 9
    @test edgecount(EtOH2) == 8

    nitro = parse(SMILES, "N#N")
    @test edgeattr(nitro, 1).order == 3

    vinylcl = parse(SMILES, "C=CCl")
    @test edgeattr(vinylcl, 1).order == 2
    @test nodeattr(vinylcl, 3).symbol == :Cl
end

@testset "branch" begin
    tBuOH = parse(SMILES, "CC(C)(C)O")
    @test nodeattr(tBuOH, 5).symbol == :O
    @test degree(tBuOH, 2) == 4

    TFA = parse(SMILES, "FC(F)(F)C(=O)O")
    @test edgeattr(TFA, 5).order == 2
    @test degree(TFA, 2) == 4

    EDTA = parse(SMILES, "C(CN(CC(=O)O)CC(=O)O)N(CC(=O)O)CC(=O)O")
    @test degree(EDTA, 5) == 3
    @test degree(EDTA, 12) == 3
    @test degree(EDTA, 14) == 3
end

@testset "ring" begin
    cychexe = parse(SMILES, "C1C=CCCC1")
    @test degree(cychexe, 1) == 2

    cychexe2 = parse(SMILES, "C=1CCCCC1")
    @test edgeattr(cychexe2, 6).order == 2

    cychexe3 = parse(SMILES, "C1CCCCC=1")
    @test edgeattr(cychexe3, 6).order == 2

    cychexe4 = parse(SMILES, "C=1CCCCC=1")
    @test edgeattr(cychexe4, 6).order == 2

    bicyclohexyl = parse(SMILES, "C1CCCCC1C1CCCCC1")
    @test edgeattr(bicyclohexyl, 6).order == 1
    @test edgeattr(bicyclohexyl, 13).order == 1

    bicyclohexyl2 = parse(SMILES, "C1CCCCC1C2CCCCC2")
    @test edgeattr(bicyclohexyl2, 6).order == 1
    @test edgeattr(bicyclohexyl2, 13).order == 1

    naphthalene = parse(SMILES, "C=1C=CC=C2C1C=CC=C2")
    @test edgeattr(naphthalene, 6).order == 2
    @test edgeattr(naphthalene, 11).order == 1

    naphthalene2 = parse(SMILES, "C1=CC=C2C(=C1)C=CC=C2")
    @test edgeattr(naphthalene2, 6).order == 1
    @test edgeattr(naphthalene2, 11).order == 1

    tetrahedrane = parse(SMILES, "C12C3C1C23")
    @test edgecount(tetrahedrane) == 6

    thiamin = parse(SMILES, "OCCc1c(C)[n+](=cs1)Cc2cnc(C)nc(N)2")
    @test edgecount(thiamin) == 19
end

@testset "component" begin
    NaCl = parse(SMILES, "[Na+].[Cl-]")
    @test nodeattr(NaCl, 1).symbol == :Na
    @test nodeattr(NaCl, 1).charge == 1
    @test nodeattr(NaCl, 2).symbol == :Cl
    @test nodeattr(NaCl, 2).charge == -1

    ferrous = parse(SMILES, "[Fe++]")
    @test nodeattr(ferrous, 1).symbol == :Fe
    @test nodeattr(ferrous, 1).charge == 2

    ferrous2 = parse(SMILES, "[Fe+2]")
    @test nodeattr(ferrous2, 1).symbol == :Fe
    @test nodeattr(ferrous2, 1).charge == 2

    CuSO4 = parse(SMILES, "[Cu+2].[O-]S(=O)(=O)[O-]")
    @test nodeattr(CuSO4, 1).charge == 2
    @test nodeattr(CuSO4, 6).charge == -1
    @test degree(CuSO4, 3) == 4

    PhONa = parse(SMILES, "c1cc([O-].[Na+])ccc1")
    @test degree(PhONa, 3) == 3
    @test degree(PhONa, 5) == 0
end

@testset "stereo" begin
    LAla = parse(SMILES, "N[C@H](C)C(=O)O")
    @test degree(LAla, 2) == 4
    @test nodeattr(LAla, 2).stereo == :anticlockwise
    DAla = parse(SMILES, "N[C@@H](C)C(=O)O")
    @test degree(DAla, 2) == 4
    @test nodeattr(DAla, 2).stereo == :clockwise
end

@testset "aromatic" begin
    benzene = parse(SMILES, "c1ccccc1")
    @test nodeattr(benzene, 1).isaromatic
    @test edgeattr(benzene, 6).order == 1
    pyrrole = parse(SMILES, "[nH]1cccc1")  # pyrrole H cannot be omitted
    @test nodeattr(pyrrole, 1).isaromatic
    @test degree(pyrrole, 1) == 3
end

end # smarts.smiles


text:
function Base.show(io::IO, csinfo::CSInfo)
    println("Info: ")
    for name in fieldnames(CSInfo)
        println(io, "$name = $(getfield(csinfo, name))")
    end
end

function compress_var_string(variable::CS.Variable)
    if CS.isfixed(variable)
        return string(CS.value(variable))
    end

    sorted_vals = sort(CS.values(variable))
    if sorted_vals[1]+length(CS.values(variable))-1 == sorted_vals[end]
        return "$(sorted_vals[1]):$(sorted_vals[end])"
    end
    return string(sort(CS.values(variable)))
end

function get_str_repr(variables::Array{CS.Variable})
    if length(size(variables)) == 1
        output = ""
        for i=1:length(variables)
            if !CS.isfixed(variables[i])
                output *= "$(compress_var_string(variables[i])), "
            else
                output *= "$(CS.value(variables[i])), "
            end
        end
        return [output[1:end-2]]
    elseif length(size(variables)) == 2
        max_length = 1
        for j=1:size(variables)[2]
            for i=1:size(variables)[1]
                if !CS.isfixed(variables[i,j])
                    len = length(compress_var_string(variables[i,j]))
                    if len > max_length
                        max_length = len
                    end
                else
                    len = length(string(CS.value(variables[i,j])))
                    if len > max_length
                        max_length = len
                    end
                end
            end
        end
        lines = String[]
        for j=1:size(variables)[2]
            line = ""
            for i=1:size(variables)[1]
                pstr = ""
                if !CS.isfixed(variables[i,j])
                    pstr = compress_var_string(variables[i,j])
                else
                    pstr = string(CS.value(variables[i,j]))
                end
                space_left  = floor(Int, (max_length-length(pstr))/2)
                space_right = ceil(Int, (max_length-length(pstr))/2)
                line *= repeat(" ", space_left)*pstr*repeat(" ", space_right)*" "
            end
            push!(lines, line)
        end
        return lines
    else
        @warn "Currently not supported to print more than 2 dimensions. Maybe file an issue with your problem and desired output"
    end
end

function Base.show(io::IO, variables::Array{CS.Variable})
    lines = get_str_repr(variables)
    for line in lines
        println(line)
    end
end


text:
# Feature specifications:

# Individual dense features are represented with strings like "s1w".
typealias DFeature AbstractString

# DFvec is a String vector specifying a set of dense features.
typealias DFvec{T<:DFeature} Vector{T}

# Individual sparse features are specified using a vector of feature
# names indicating a conjunction, e.g. ["s0w", "s0p"]
typealias SFeature{T<:AbstractString} Vector{T}

# SFVec is an vector of SFeatures specifying a set of sparse features.
# e.g. [["sw"],["nw","np"]]
typealias SFvec{T<:SFeature} Vector{T}

# Fvec is a union of DFvec and SFvec, which specify dense and sparse
# features respectively.  The user controls what kind of features are
# generated by using the appropriate type of feature vector.
typealias Fvec Union{DFvec,SFvec}

"""

    features(p::Parser, s::Sentence, feats, [x, idx])

Return a feature vector given a [`Parser`](@ref) state, a
[`Sentence`](@ref), and a feature specification vector.  The last two
arguments (optional) can provide a preallocated output array and an
offset position in that array.

Primitive features are represented by strings of the following form:

    [sn]\d?([hlr]\d?)*[wvcpdLabAB]

The meaning of each character is described below. The letter i
represents an optional single digit integer. A default value is used
if i is not specified.

|char|meaning|
|:---|:------|
|si|i'th stack word, default i=0 means top|
|ni|i'th buffer word, default i=0 means first|
|hi|i'th degree head, default i=1 means direct head|
|li|i'th leftmost child, default i=1 means the leftmost child|
|ri|i'th rightmost child, default i=1 means the rightmost child|
|w|word|
|v|word vector (first half of wvec)|
|c|context vector (second half of wvec)|
|p|postag|
|d|distance to the right using encoding: 1,2,3,4,5-10,10+|
|L|dependency label (0 is ROOT or NONE)|
|a|number of left children.|
|b|number of right children.|
|A|set of left dependency labels|
|B|set of right dependency labels|

For example "sp" or "s0p" indicate the postag of the top stack word,
"n0lL" means the dependency label of the leftmost child of the first
buffer word, "s1d" is the s0s1 distance, "s0d" is the s0n0 distance.

The returned feature vector can be dense (`Vector{Float}`) or sparse
(`Vector{Int}` indicating non-zero positions in a sparse binary
vector).  The choice is marked by the type of the `feats` argument:
for dense output it should be a `Vector{String}`, for sparse output it
should be a `Vector{Vector{String}}`:

    ["s0v","n0lp","n1p"]  # dense feature spec
    [["sw"],["sw","sp"]]  # sparse feature spec

This allows sparse feature specifications to contain elements with
multiple primitive features indicating conjunctions. Some other
notable differences between dense and sparse features are:

* Dense vectors do not support "word", sparse vectors do not support "word vector" or "context vector".
* Sparse vectors represent number of children exactly, dense vectors use the encoding 0,1,...,8,9+.
* For dense features, `x` is a float matrix and `idx=1` is the column number to fill. For sparse features `x` is an int vector and `idx=0` is an offset.

"""
function features(p::Parser, s::Sentence, feats::DFvec,
                  x::AbstractMatrix=Array(wtype(s),flen(p,s,feats),1),
                  xcol::Integer=1)
    if xcol > size(x,2); error("xcol > size(x,2)"); end
    wrows = wdim(s)             #11
    xrows = size(x,1)
    xtype = eltype(x)
    x1 = one(xtype)
    x[:,xcol] = zero(xtype)     #173
    nx = 0                      # last entry in x
    nv = wrows # >> 1             # size of word/context vector, assumes word vec and context vec concatenated
    nd = p.ndeps
    np = length(s.vocab.postags)
    nw = p.nword
    ldep,rdep,lset,rset = getdeps(p) #3775
    for f in feats                   #48
        a = getanchor(f,p,ldep,rdep) #1402
        fn = f[end]                  #642
        if fn == 'v'                 #234
            if a>0                   #91
                copy!(x, (xcol-1)*xrows+nx+1, s.vocab.wvecs, (s.word[a]-1)*wrows+1, nv) #331
            end; nx += nv       #5
        elseif fn == 'c'        #200
            error("Context vectors not implemented.")
            # if a>0              #122
            #     copy!(x, (xcol-1)*xrows+nx+1, s.vocab.wvec, (a-1)*wrows+nv+1, nv) #492
            # end; nx += nv       #6
        elseif fn == 'p'        #114
            if a>0              #159
                if s.postag[a] > np; error("postag out of bound"); end #210
                x[nx+s.postag[a], xcol] = x1 #839
            end; nx += np                    #5
        elseif fn == 'd'
            if a>0
                d = getrdist(f,p,a)
                if d>0; x[nx+(d>10?6:d>5?5:d), xcol] = x1; end
            end; nx += 6
        elseif fn == 'L'
            if a>0
                r = p.deprel[a]
                if r > nd; error("deprel out of bound"); end
                x[nx+1+r, xcol] = x1 # first bit for deprel=0 (ROOT)
            end; nx += (nd+1)
        elseif fn == 'a'
            if a>0
                if isassigned(ldep,a)
                    lcnt=length(ldep[a])
                    if lcnt > 9; lcnt = 9; end
                else
                    lcnt = 0
                end
                x[nx+1+lcnt, xcol] = x1 # 0-9
            end; nx += 10
        elseif fn == 'b'
            if a>0
                if isassigned(ldep,a)
                    rcnt=length(ldep[a])
                    if rcnt > 9; rcnt = 9; end
                else
                    rcnt = 0
                end
                x[nx+1+rcnt, xcol] = x1 # 0-9
            end; nx += 10
        elseif fn == 'A'
            if a>0 && isassigned(lset,a)
                copy!(x, (xcol-1)*xrows+nx+1, Array{wtype(s)}(lset[a]), 1, nd) #
            end; nx += nd
        elseif fn == 'B'
            if a>0 && isassigned(rset,a)
                copy!(x, (xcol-1)*xrows+nx+1, Array{wtype(s)}(rset[a]), 1, nd) #
            end; nx += nd
        elseif fn == 'w'
            error("Dense features do not support 'w'")
        else
            error("Unknown feature $(fn)") # 3
        end
    end
    nx == xrows || error("Bad feature vector length $nx != $xrows")
    return x
end

# Sparse feature extractor:

# TODO: change default idx=1.

function features(p::Parser, s::Sentence, feats::SFvec,
                  x::AbstractMatrix=Array(SFtype, length(feats), 1), xcol=1)
    if xcol > size(x,2); error("xcol > size(x,2)"); end
    deps = getdeps(p)
    SFhash = s.vocab.fdict
    @inbounds for i = 1:length(feats)
        f = feats[i]
        v = Array(Any, length(f))                             # TODO: get rid of alloc here?
        for j=1:length(f)
            v[j] = features1(p,s,f[j],deps...) #1023
        end
        x[i,xcol] = get!(SFhash, (f,v), 1+length(SFhash)) #7782 TODO: need a better hash function here
    end
    return x
end

# Here is where the actual feature lookup happens.  Similar to the
# dense lookup.  But does not have to convert everything to a number,
# can return strings etc.  Returns 'nothing' if target word does not
# exist or the feature is not available.

function features1(p::Parser, s::Sentence, f::String, ldep, rdep, lset, rset)
    a = getanchor(f,p,ldep,rdep)
    if a == 0; return nothing; end
    fn = f[end]
    if fn == 'w'; s.word[a]
    elseif fn == 'p'; s.postag[a]
    elseif fn == 'L'; p.deprel[a]
    elseif fn == 'a'; if isassigned(ldep,a); length(ldep[a]); else; 0; end
    elseif fn == 'b'; if isassigned(rdep,a); length(rdep[a]); else; 0; end
    elseif fn == 'A'; if isassigned(lset,a); lset[a]; else; nothing; end
    elseif fn == 'B'; if isassigned(rset,a); rset[a]; else; nothing; end
    elseif fn == 'd'; (d=getrdist(f,p,a); if d>10; 6; elseif d>5; 5; else; d; end)
    elseif fn == 'v'; error("Sparse features do not support 'v'")
    elseif fn == 'c'; error("Sparse features do not support 'c'")
    else error("Unknown feature letter $fn")
    end
end

function flen(p::Parser, s::Sentence, feats::SFvec)
    length(feats) # one integer per feature
end

function flen(p::Parser, s::Sentence, feats::DFvec)
    nx = 0
    nw = wdim(s)  # >> 1
    nd = p.ndeps
    np = length(s.vocab.postags)
    for f in feats
        nx += flen1(f[end], nw, nd, np) # 1129
    end
    return nx
end

function flen1(c::Char, nw::Int, nd::Int, np::Int)
    if c == 'v'; nw
    elseif c == 'c'; nw
    elseif c == 'p'; np
    elseif c == 'd'; 6
    elseif c == 'L'; (nd+1) # ROOT is not included in nd, TODO: fix this
    elseif c == 'a'; 10
    elseif c == 'b'; 10
    elseif c == 'A'; nd
    elseif c == 'B'; nd
    else error("Unknown feature character $c")
    end
end

# Utility functions to calculate the size of the feature matrix

xsize(p::Parser, s::Sentence, f::Fvec)=(flen(p,s,f),nmoves(p,s))
xsize(p::Parser, c::Corpus, f::Fvec)=(flen(p,c[1],f),nmoves(p,c))
xsize{T<:Parser}(p::Vector{T}, c::Corpus, f::Fvec)=xsize(p[1],c,f)
ysize(p::Parser, s::Sentence)=(nmoves(p),nmoves(p,s))
ysize(p::Parser, c::Corpus)=(nmoves(p),nmoves(p,c))
ysize{T<:Parser}(p::Vector{T}, c::Corpus)=ysize(p[1],c)

function getdeps{T<:Parser}(p::T)
    nw = p.nword
    nd = p.ndeps
    dr = p.deprel
    ldep = Array(Any,nw)
    rdep = Array(Any,nw)
    lset = Array(Any,nw)
    rset = Array(Any,nw)
    @inbounds for d=1:nw
        h=Int(p.head[d])
        if h==0
            continue
        elseif d<h
            if !isassigned(ldep,h)
                ldep[h]=[d]
                lset[h]=zeros(UInt8,nd) # falses(nd) is slower, zeros(WVtype,nd) same speed more space
                lset[h][dr[d]]=1
            else
                push!(ldep[h],d)
                lset[h][dr[d]]=1
            end
        elseif d>h
            if !isassigned(rdep,h)
                rdep[h]=[d]
                rset[h]=zeros(UInt8,nd)
                rset[h][dr[d]]=1
            else
                unshift!(rdep[h],d)
                rset[h][dr[d]]=1
            end
        else
            error("h==d")
        end
    end
    return (ldep,rdep,lset,rset)
end

function getanchor(f::String, p::Parser, ldep, rdep)
    f1 = f[1]; f2 = f[2]; flen = length(f)
    if isdigit(f2)
        i = f2 - '0'            # index
        n = 3                   # next character
    else
        i = 0
        n = 2
    end
    a = 0                       # target word
    if f1 == 's'
        if (p.sptr - i >= 1) # 25
            a = Int(p.stack[p.sptr - i])             # 456
        end
    elseif f1 == 'n'
        if p.wptr + i <= p.nword
            a = p.wptr + i
        end
    else
        error("feature string should start with [sn]")
    end
    if a==0; return 0; end
    while n < flen
        f1 = f[n]; f2 = f[n+1]
        if isdigit(f2)
            i = f2 - '0'
            n = n+2
        else
            i = 1
            n = n+1
        end
        if i <= 0
            error("hlr indexing is one based") # 3 [lrh] is one based, [sn] was zero based
        end
        if f1 == 'l'                              # 2
            if a > p.wptr; error("buffer words other than n0 do not have ldeps"); end # 252
            if isassigned(ldep,a) && i <= length(ldep[a])
                a = Int(ldep[a][i])
            else
                return 0
            end
        elseif f1 == 'r'
            if a >= p.wptr; error("buffer words do not have rdeps"); end
            if isassigned(rdep,a) && i <= length(rdep[a])
                a = Int(rdep[a][i])
            else
                return 0
            end
        elseif f1 == 'h'
            for j=1:i       # 5
                a = Int(p.head[a]) # 147
                if a == 0
                    return 0
                end
            end
        else
            break
        end
    end
    if n != length(f); error("n!=length(f)"); end
    return a
end

function getrdist(f::AbstractString, p::Parser, a::Integer)
    if f[1]=='s'
        if isdigit(f[2])
            i=f[2]-'0'
        else
            i=0
        end
        if i>0
            return p.stack[p.sptr - i + 1] - a
        elseif p.wptr <= p.nword
            return p.wptr - a
        else
            return 0
        end
    else
        return 0 # dist only defined for stack words
    end
end


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

"""
    activate_anomaly_detector(anomaly_detector_arn)
    activate_anomaly_detector(anomaly_detector_arn, params::Dict{String,<:Any})

Activates an anomaly detector.

# Arguments
- `anomaly_detector_arn`: The ARN of the anomaly detector.

"""
activate_anomaly_detector(AnomalyDetectorArn; aws_config::AbstractAWSConfig=global_aws_config()) = lookoutmetrics("POST", "/ActivateAnomalyDetector", Dict{String, Any}("AnomalyDetectorArn"=>AnomalyDetectorArn); aws_config=aws_config)
activate_anomaly_detector(AnomalyDetectorArn, params::AbstractDict{String}; aws_config::AbstractAWSConfig=global_aws_config()) = lookoutmetrics("POST", "/ActivateAnomalyDetector", Dict{String, Any}(mergewith(_merge, Dict{String, Any}("AnomalyDetectorArn"=>AnomalyDetectorArn), params)); aws_config=aws_config)

"""
    back_test_anomaly_detector(anomaly_detector_arn)
    back_test_anomaly_detector(anomaly_detector_arn, params::Dict{String,<:Any})

Runs a backtest for anomaly detection for the specified resource.

# Arguments
- `anomaly_detector_arn`: The Amazon Resource Name (ARN) of the anomaly detector.

"""
back_test_anomaly_detector(AnomalyDetectorArn; aws_config::AbstractAWSConfig=global_aws_config()) = lookoutmetrics("POST", "/BackTestAnomalyDetector", Dict{String, Any}("AnomalyDetectorArn"=>AnomalyDetectorArn); aws_config=aws_config)
back_test_anomaly_detector(AnomalyDetectorArn, params::AbstractDict{String}; aws_config::AbstractAWSConfig=global_aws_config()) = lookoutmetrics("POST", "/BackTestAnomalyDetector", Dict{String, Any}(mergewith(_merge, Dict{String, Any}("AnomalyDetectorArn"=>AnomalyDetectorArn), params)); aws_config=aws_config)

"""
    create_alert(action, alert_name, alert_sensitivity_threshold, anomaly_detector_arn)
    create_alert(action, alert_name, alert_sensitivity_threshold, anomaly_detector_arn, params::Dict{String,<:Any})

Creates an alert for an anomaly detector.

# Arguments
- `action`: Action that will be triggered when there is an alert.
- `alert_name`: The name of the alert.
- `alert_sensitivity_threshold`: An integer from 0 to 100 specifying the alert sensitivity
  threshold.
- `anomaly_detector_arn`: The ARN of the detector to which the alert is attached.

# Optional Parameters
Optional parameters can be passed as a `params::Dict{String,<:Any}`. Valid keys are:
- `"AlertDescription"`: A description of the alert.
- `"Tags"`: A list of tags to apply to the alert.
"""
create_alert(Action, AlertName, AlertSensitivityThreshold, AnomalyDetectorArn; aws_config::AbstractAWSConfig=global_aws_config()) = lookoutmetrics("POST", "/CreateAlert", Dict{String, Any}("Action"=>Action, "AlertName"=>AlertName, "AlertSensitivityThreshold"=>AlertSensitivityThreshold, "AnomalyDetectorArn"=>AnomalyDetectorArn); aws_config=aws_config)
create_alert(Action, AlertName, AlertSensitivityThreshold, AnomalyDetectorArn, params::AbstractDict{String}; aws_config::AbstractAWSConfig=global_aws_config()) = lookoutmetrics("POST", "/CreateAlert", Dict{String, Any}(mergewith(_merge, Dict{String, Any}("Action"=>Action, "AlertName"=>AlertName, "AlertSensitivityThreshold"=>AlertSensitivityThreshold, "AnomalyDetectorArn"=>AnomalyDetectorArn), params)); aws_config=aws_config)

"""
    create_anomaly_detector(anomaly_detector_config, anomaly_detector_name)
    create_anomaly_detector(anomaly_detector_config, anomaly_detector_name, params::Dict{String,<:Any})

Creates an anomaly detector.

# Arguments
- `anomaly_detector_config`: Contains information about the configuration of the anomaly
  detector.
- `anomaly_detector_name`: The name of the detector.

# Optional Parameters
Optional parameters can be passed as a `params::Dict{String,<:Any}`. Valid keys are:
- `"AnomalyDetectorDescription"`: A description of the detector.
- `"KmsKeyArn"`: The ARN of the KMS key to use to encrypt your data.
- `"Tags"`: A list of tags to apply to the anomaly detector.
"""
create_anomaly_detector(AnomalyDetectorConfig, AnomalyDetectorName; aws_config::AbstractAWSConfig=global_aws_config()) = lookoutmetrics("POST", "/CreateAnomalyDetector", Dict{String, Any}("AnomalyDetectorConfig"=>AnomalyDetectorConfig, "AnomalyDetectorName"=>AnomalyDetectorName); aws_config=aws_config)
create_anomaly_detector(AnomalyDetectorConfig, AnomalyDetectorName, params::AbstractDict{String}; aws_config::AbstractAWSConfig=global_aws_config()) = lookoutmetrics("POST", "/CreateAnomalyDetector", Dict{String, Any}(mergewith(_merge, Dict{String, Any}("AnomalyDetectorConfig"=>AnomalyDetectorConfig, "AnomalyDetectorName"=>AnomalyDetectorName), params)); aws_config=aws_config)

"""
    create_metric_set(anomaly_detector_arn, metric_list, metric_set_name, metric_source)
    create_metric_set(anomaly_detector_arn, metric_list, metric_set_name, metric_source, params::Dict{String,<:Any})

Creates a dataset.

# Arguments
- `anomaly_detector_arn`: The ARN of the anomaly detector that will use the dataset.
- `metric_list`: A list of metrics that the dataset will contain.
- `metric_set_name`: The name of the dataset.
- `metric_source`: Contains information about how the source data should be interpreted.

# Optional Parameters
Optional parameters can be passed as a `params::Dict{String,<:Any}`. Valid keys are:
- `"DimensionList"`: A list of the fields you want to treat as dimensions.
- `"MetricSetDescription"`: A description of the dataset you are creating.
- `"MetricSetFrequency"`: The frequency with which the source data will be analyzed for
  anomalies.
- `"Offset"`: After an interval ends, the amount of time that the detector waits before
  importing data.
- `"Tags"`: A list of tags to apply to the dataset.
- `"TimestampColumn"`: Contains information about the column used for tracking time in your
  source data.
- `"Timezone"`: The time zone in which your source data was recorded.
"""
create_metric_set(AnomalyDetectorArn, MetricList, MetricSetName, MetricSource; aws_config::AbstractAWSConfig=global_aws_config()) = lookoutmetrics("POST", "/CreateMetricSet", Dict{String, Any}("AnomalyDetectorArn"=>AnomalyDetectorArn, "MetricList"=>MetricList, "MetricSetName"=>MetricSetName, "MetricSource"=>MetricSource); aws_config=aws_config)
create_metric_set(AnomalyDetectorArn, MetricList, MetricSetName, MetricSource, params::AbstractDict{String}; aws_config::AbstractAWSConfig=global_aws_config()) = lookoutmetrics("POST", "/CreateMetricSet", Dict{String, Any}(mergewith(_merge, Dict{String, Any}("AnomalyDetectorArn"=>AnomalyDetectorArn, "MetricList"=>MetricList, "MetricSetName"=>MetricSetName, "MetricSource"=>MetricSource), params)); aws_config=aws_config)

"""
    delete_alert(alert_arn)
    delete_alert(alert_arn, params::Dict{String,<:Any})

Deletes an alert.

# Arguments
- `alert_arn`: The ARN of the alert to delete.

"""
delete_alert(AlertArn; aws_config::AbstractAWSConfig=global_aws_config()) = lookoutmetrics("POST", "/DeleteAlert", Dict{String, Any}("AlertArn"=>AlertArn); aws_config=aws_config)
delete_alert(AlertArn, params::AbstractDict{String}; aws_config::AbstractAWSConfig=global_aws_config()) = lookoutmetrics("POST", "/DeleteAlert", Dict{String, Any}(mergewith(_merge, Dict{String, Any}("AlertArn"=>AlertArn), params)); aws_config=aws_config)

"""
    delete_anomaly_detector(anomaly_detector_arn)
    delete_anomaly_detector(anomaly_detector_arn, params::Dict{String,<:Any})

Deletes a detector. Deleting an anomaly detector will delete all of its corresponding
resources including any configured datasets and alerts.

# Arguments
- `anomaly_detector_arn`: The ARN of the detector to delete.

"""
delete_anomaly_detector(AnomalyDetectorArn; aws_config::AbstractAWSConfig=global_aws_config()) = lookoutmetrics("POST", "/DeleteAnomalyDetector", Dict{String, Any}("AnomalyDetectorArn"=>AnomalyDetectorArn); aws_config=aws_config)
delete_anomaly_detector(AnomalyDetectorArn, params::AbstractDict{String}; aws_config::AbstractAWSConfig=global_aws_config()) = lookoutmetrics("POST", "/DeleteAnomalyDetector", Dict{String, Any}(mergewith(_merge, Dict{String, Any}("AnomalyDetectorArn"=>AnomalyDetectorArn), params)); aws_config=aws_config)

"""
    describe_alert(alert_arn)
    describe_alert(alert_arn, params::Dict{String,<:Any})

Describes an alert. Amazon Lookout for Metrics API actions are eventually consistent. If
you do a read operation on a resource immediately after creating or modifying it, use
retries to allow time for the write operation to complete.

# Arguments
- `alert_arn`: The ARN of the alert to describe.

"""
describe_alert(AlertArn; aws_config::AbstractAWSConfig=global_aws_config()) = lookoutmetrics("POST", "/DescribeAlert", Dict{String, Any}("AlertArn"=>AlertArn); aws_config=aws_config)
describe_alert(AlertArn, params::AbstractDict{String}; aws_config::AbstractAWSConfig=global_aws_config()) = lookoutmetrics("POST", "/DescribeAlert", Dict{String, Any}(mergewith(_merge, Dict{String, Any}("AlertArn"=>AlertArn), params)); aws_config=aws_config)

"""
    describe_anomaly_detection_executions(anomaly_detector_arn)
    describe_anomaly_detection_executions(anomaly_detector_arn, params::Dict{String,<:Any})

Returns information about the status of the specified anomaly detection jobs.

# Arguments
- `anomaly_detector_arn`: The Amazon Resource Name (ARN) of the anomaly detector.

# Optional Parameters
Optional parameters can be passed as a `params::Dict{String,<:Any}`. Valid keys are:
- `"MaxResults"`: The number of items to return in the response.
- `"NextToken"`: Specify the pagination token that's returned by a previous request to
  retrieve the next page of results.
- `"Timestamp"`: The timestamp of the anomaly detection job.
"""
describe_anomaly_detection_executions(AnomalyDetectorArn; aws_config::AbstractAWSConfig=global_aws_config()) = lookoutmetrics("POST", "/DescribeAnomalyDetectionExecutions", Dict{String, Any}("AnomalyDetectorArn"=>AnomalyDetectorArn); aws_config=aws_config)
describe_anomaly_detection_executions(AnomalyDetectorArn, params::AbstractDict{String}; aws_config::AbstractAWSConfig=global_aws_config()) = lookoutmetrics("POST", "/DescribeAnomalyDetectionExecutions", Dict{String, Any}(mergewith(_merge, Dict{String, Any}("AnomalyDetectorArn"=>AnomalyDetectorArn), params)); aws_config=aws_config)

"""
    describe_anomaly_detector(anomaly_detector_arn)
    describe_anomaly_detector(anomaly_detector_arn, params::Dict{String,<:Any})

Describes a detector. Amazon Lookout for Metrics API actions are eventually consistent. If
you do a read operation on a resource immediately after creating or modifying it, use
retries to allow time for the write operation to complete.

# Arguments
- `anomaly_detector_arn`: The ARN of the detector to describe.

"""
describe_anomaly_detector(AnomalyDetectorArn; aws_config::AbstractAWSConfig=global_aws_config()) = lookoutmetrics("POST", "/DescribeAnomalyDetector", Dict{String, Any}("AnomalyDetectorArn"=>AnomalyDetectorArn); aws_config=aws_config)
describe_anomaly_detector(AnomalyDetectorArn, params::AbstractDict{String}; aws_config::AbstractAWSConfig=global_aws_config()) = lookoutmetrics("POST", "/DescribeAnomalyDetector", Dict{String, Any}(mergewith(_merge, Dict{String, Any}("AnomalyDetectorArn"=>AnomalyDetectorArn), params)); aws_config=aws_config)

"""
    describe_metric_set(metric_set_arn)
    describe_metric_set(metric_set_arn, params::Dict{String,<:Any})

Describes a dataset. Amazon Lookout for Metrics API actions are eventually consistent. If
you do a read operation on a resource immediately after creating or modifying it, use
retries to allow time for the write operation to complete.

# Arguments
- `metric_set_arn`: The ARN of the dataset.

"""
describe_metric_set(MetricSetArn; aws_config::AbstractAWSConfig=global_aws_config()) = lookoutmetrics("POST", "/DescribeMetricSet", Dict{String, Any}("MetricSetArn"=>MetricSetArn); aws_config=aws_config)
describe_metric_set(MetricSetArn, params::AbstractDict{String}; aws_config::AbstractAWSConfig=global_aws_config()) = lookoutmetrics("POST", "/DescribeMetricSet", Dict{String, Any}(mergewith(_merge, Dict{String, Any}("MetricSetArn"=>MetricSetArn), params)); aws_config=aws_config)

"""
    get_anomaly_group(anomaly_detector_arn, anomaly_group_id)
    get_anomaly_group(anomaly_detector_arn, anomaly_group_id, params::Dict{String,<:Any})

Returns details about a group of anomalous metrics.

# Arguments
- `anomaly_detector_arn`: The Amazon Resource Name (ARN) of the anomaly detector.
- `anomaly_group_id`: The ID of the anomaly group.

"""
get_anomaly_group(AnomalyDetectorArn, AnomalyGroupId; aws_config::AbstractAWSConfig=global_aws_config()) = lookoutmetrics("POST", "/GetAnomalyGroup", Dict{String, Any}("AnomalyDetectorArn"=>AnomalyDetectorArn, "AnomalyGroupId"=>AnomalyGroupId); aws_config=aws_config)
get_anomaly_group(AnomalyDetectorArn, AnomalyGroupId, params::AbstractDict{String}; aws_config::AbstractAWSConfig=global_aws_config()) = lookoutmetrics("POST", "/GetAnomalyGroup", Dict{String, Any}(mergewith(_merge, Dict{String, Any}("AnomalyDetectorArn"=>AnomalyDetectorArn, "AnomalyGroupId"=>AnomalyGroupId), params)); aws_config=aws_config)

"""
    get_feedback(anomaly_detector_arn, anomaly_group_time_series_feedback)
    get_feedback(anomaly_detector_arn, anomaly_group_time_series_feedback, params::Dict{String,<:Any})

Get feedback for an anomaly group.

# Arguments
- `anomaly_detector_arn`: The Amazon Resource Name (ARN) of the anomaly detector.
- `anomaly_group_time_series_feedback`: The anomalous metric and group ID.

# Optional Parameters
Optional parameters can be passed as a `params::Dict{String,<:Any}`. Valid keys are:
- `"MaxResults"`: The maximum number of results to return.
- `"NextToken"`: Specify the pagination token that's returned by a previous request to
  retrieve the next page of results.
"""
get_feedback(AnomalyDetectorArn, AnomalyGroupTimeSeriesFeedback; aws_config::AbstractAWSConfig=global_aws_config()) = lookoutmetrics("POST", "/GetFeedback", Dict{String, Any}("AnomalyDetectorArn"=>AnomalyDetectorArn, "AnomalyGroupTimeSeriesFeedback"=>AnomalyGroupTimeSeriesFeedback); aws_config=aws_config)
get_feedback(AnomalyDetectorArn, AnomalyGroupTimeSeriesFeedback, params::AbstractDict{String}; aws_config::AbstractAWSConfig=global_aws_config()) = lookoutmetrics("POST", "/GetFeedback", Dict{String, Any}(mergewith(_merge, Dict{String, Any}("AnomalyDetectorArn"=>AnomalyDetectorArn, "AnomalyGroupTimeSeriesFeedback"=>AnomalyGroupTimeSeriesFeedback), params)); aws_config=aws_config)

"""
    get_sample_data()
    get_sample_data(params::Dict{String,<:Any})

Returns a selection of sample records from an Amazon S3 datasource.

# Optional Parameters
Optional parameters can be passed as a `params::Dict{String,<:Any}`. Valid keys are:
- `"S3SourceConfig"`: A datasource bucket in Amazon S3.
"""
get_sample_data(; aws_config::AbstractAWSConfig=global_aws_config()) = lookoutmetrics("POST", "/GetSampleData"; aws_config=aws_config)
get_sample_data(params::AbstractDict{String}; aws_config::AbstractAWSConfig=global_aws_config()) = lookoutmetrics("POST", "/GetSampleData", params; aws_config=aws_config)

"""
    list_alerts()
    list_alerts(params::Dict{String,<:Any})

Lists the alerts attached to a detector. Amazon Lookout for Metrics API actions are
eventually consistent. If you do a read operation on a resource immediately after creating
or modifying it, use retries to allow time for the write operation to complete.

# Optional Parameters
Optional parameters can be passed as a `params::Dict{String,<:Any}`. Valid keys are:
- `"AnomalyDetectorArn"`: The ARN of the alert's detector.
- `"MaxResults"`: The maximum number of results that will be displayed by the request.
- `"NextToken"`: If the result of the previous request is truncated, the response includes
  a NextToken. To retrieve the next set of results, use the token in the next request. Tokens
  expire after 24 hours.
"""
list_alerts(; aws_config::AbstractAWSConfig=global_aws_config()) = lookoutmetrics("POST", "/ListAlerts"; aws_config=aws_config)
list_alerts(params::AbstractDict{String}; aws_config::AbstractAWSConfig=global_aws_config()) = lookoutmetrics("POST", "/ListAlerts", params; aws_config=aws_config)

"""
    list_anomaly_detectors()
    list_anomaly_detectors(params::Dict{String,<:Any})

Lists the detectors in the current AWS Region. Amazon Lookout for Metrics API actions are
eventually consistent. If you do a read operation on a resource immediately after creating
or modifying it, use retries to allow time for the write operation to complete.

# Optional Parameters
Optional parameters can be passed as a `params::Dict{String,<:Any}`. Valid keys are:
- `"MaxResults"`: The maximum number of results to return.
- `"NextToken"`: If the result of the previous request was truncated, the response includes
  a NextToken. To retrieve the next set of results, use the token in the next request. Tokens
  expire after 24 hours.
"""
list_anomaly_detectors(; aws_config::AbstractAWSConfig=global_aws_config()) = lookoutmetrics("POST", "/ListAnomalyDetectors"; aws_config=aws_config)
list_anomaly_detectors(params::AbstractDict{String}; aws_config::AbstractAWSConfig=global_aws_config()) = lookoutmetrics("POST", "/ListAnomalyDetectors", params; aws_config=aws_config)

"""
    list_anomaly_group_summaries(anomaly_detector_arn, sensitivity_threshold)
    list_anomaly_group_summaries(anomaly_detector_arn, sensitivity_threshold, params::Dict{String,<:Any})

Returns a list of anomaly groups.

# Arguments
- `anomaly_detector_arn`: The Amazon Resource Name (ARN) of the anomaly detector.
- `sensitivity_threshold`: The minimum severity score for inclusion in the output.

# Optional Parameters
Optional parameters can be passed as a `params::Dict{String,<:Any}`. Valid keys are:
- `"MaxResults"`: The maximum number of results to return.
- `"NextToken"`: Specify the pagination token that's returned by a previous request to
  retrieve the next page of results.
"""
list_anomaly_group_summaries(AnomalyDetectorArn, SensitivityThreshold; aws_config::AbstractAWSConfig=global_aws_config()) = lookoutmetrics("POST", "/ListAnomalyGroupSummaries", Dict{String, Any}("AnomalyDetectorArn"=>AnomalyDetectorArn, "SensitivityThreshold"=>SensitivityThreshold); aws_config=aws_config)
list_anomaly_group_summaries(AnomalyDetectorArn, SensitivityThreshold, params::AbstractDict{String}; aws_config::AbstractAWSConfig=global_aws_config()) = lookoutmetrics("POST", "/ListAnomalyGroupSummaries", Dict{String, Any}(mergewith(_merge, Dict{String, Any}("AnomalyDetectorArn"=>AnomalyDetectorArn, "SensitivityThreshold"=>SensitivityThreshold), params)); aws_config=aws_config)

"""
    list_anomaly_group_time_series(anomaly_detector_arn, anomaly_group_id, metric_name)
    list_anomaly_group_time_series(anomaly_detector_arn, anomaly_group_id, metric_name, params::Dict{String,<:Any})

Gets a list of anomalous metrics for a measure in an anomaly group.

# Arguments
- `anomaly_detector_arn`: The Amazon Resource Name (ARN) of the anomaly detector.
- `anomaly_group_id`: The ID of the anomaly group.
- `metric_name`: The name of the measure field.

# Optional Parameters
Optional parameters can be passed as a `params::Dict{String,<:Any}`. Valid keys are:
- `"MaxResults"`: The maximum number of results to return.
- `"NextToken"`: Specify the pagination token that's returned by a previous request to
  retrieve the next page of results.
"""
list_anomaly_group_time_series(AnomalyDetectorArn, AnomalyGroupId, MetricName; aws_config::AbstractAWSConfig=global_aws_config()) = lookoutmetrics("POST", "/ListAnomalyGroupTimeSeries", Dict{String, Any}("AnomalyDetectorArn"=>AnomalyDetectorArn, "AnomalyGroupId"=>AnomalyGroupId, "MetricName"=>MetricName); aws_config=aws_config)
list_anomaly_group_time_series(AnomalyDetectorArn, AnomalyGroupId, MetricName, params::AbstractDict{String}; aws_config::AbstractAWSConfig=global_aws_config()) = lookoutmetrics("POST", "/ListAnomalyGroupTimeSeries", Dict{String, Any}(mergewith(_merge, Dict{String, Any}("AnomalyDetectorArn"=>AnomalyDetectorArn, "AnomalyGroupId"=>AnomalyGroupId, "MetricName"=>MetricName), params)); aws_config=aws_config)

"""
    list_metric_sets()
    list_metric_sets(params::Dict{String,<:Any})

Lists the datasets in the current AWS Region. Amazon Lookout for Metrics API actions are
eventually consistent. If you do a read operation on a resource immediately after creating
or modifying it, use retries to allow time for the write operation to complete.

# Optional Parameters
Optional parameters can be passed as a `params::Dict{String,<:Any}`. Valid keys are:
- `"AnomalyDetectorArn"`: The ARN of the anomaly detector containing the metrics sets to
  list.
- `"MaxResults"`: The maximum number of results to return.
- `"NextToken"`: If the result of the previous request was truncated, the response includes
  a NextToken. To retrieve the next set of results, use the token in the next request. Tokens
  expire after 24 hours.
"""
list_metric_sets(; aws_config::AbstractAWSConfig=global_aws_config()) = lookoutmetrics("POST", "/ListMetricSets"; aws_config=aws_config)
list_metric_sets(params::AbstractDict{String}; aws_config::AbstractAWSConfig=global_aws_config()) = lookoutmetrics("POST", "/ListMetricSets", params; aws_config=aws_config)

"""
    list_tags_for_resource(resource_arn)
    list_tags_for_resource(resource_arn, params::Dict{String,<:Any})

Gets a list of tags for a detector, dataset, or alert.

# Arguments
- `resource_arn`: The resource's Amazon Resource Name (ARN).

"""
list_tags_for_resource(resourceArn; aws_config::AbstractAWSConfig=global_aws_config()) = lookoutmetrics("GET", "/tags/$(resourceArn)"; aws_config=aws_config)
list_tags_for_resource(resourceArn, params::AbstractDict{String}; aws_config::AbstractAWSConfig=global_aws_config()) = lookoutmetrics("GET", "/tags/$(resourceArn)", params; aws_config=aws_config)

"""
    put_feedback(anomaly_detector_arn, anomaly_group_time_series_feedback)
    put_feedback(anomaly_detector_arn, anomaly_group_time_series_feedback, params::Dict{String,<:Any})

Add feedback for an anomalous metric.

# Arguments
- `anomaly_detector_arn`: The Amazon Resource Name (ARN) of the anomaly detector.
- `anomaly_group_time_series_feedback`: Feedback for an anomalous metric.

"""
put_feedback(AnomalyDetectorArn, AnomalyGroupTimeSeriesFeedback; aws_config::AbstractAWSConfig=global_aws_config()) = lookoutmetrics("POST", "/PutFeedback", Dict{String, Any}("AnomalyDetectorArn"=>AnomalyDetectorArn, "AnomalyGroupTimeSeriesFeedback"=>AnomalyGroupTimeSeriesFeedback); aws_config=aws_config)
put_feedback(AnomalyDetectorArn, AnomalyGroupTimeSeriesFeedback, params::AbstractDict{String}; aws_config::AbstractAWSConfig=global_aws_config()) = lookoutmetrics("POST", "/PutFeedback", Dict{String, Any}(mergewith(_merge, Dict{String, Any}("AnomalyDetectorArn"=>AnomalyDetectorArn, "AnomalyGroupTimeSeriesFeedback"=>AnomalyGroupTimeSeriesFeedback), params)); aws_config=aws_config)

"""
    tag_resource(resource_arn, tags)
    tag_resource(resource_arn, tags, params::Dict{String,<:Any})

Adds tags to a detector, dataset, or alert.

# Arguments
- `resource_arn`: The resource's Amazon Resource Name (ARN).
- `tags`: Tags to apply to the resource. Tag keys and values can contain letters, numbers,
  spaces, and the following symbols: _.:/=+@-

"""
tag_resource(resourceArn, tags; aws_config::AbstractAWSConfig=global_aws_config()) = lookoutmetrics("POST", "/tags/$(resourceArn)", Dict{String, Any}("tags"=>tags); aws_config=aws_config)
tag_resource(resourceArn, tags, params::AbstractDict{String}; aws_config::AbstractAWSConfig=global_aws_config()) = lookoutmetrics("POST", "/tags/$(resourceArn)", Dict{String, Any}(mergewith(_merge, Dict{String, Any}("tags"=>tags), params)); aws_config=aws_config)

"""
    untag_resource(resource_arn, tag_keys)
    untag_resource(resource_arn, tag_keys, params::Dict{String,<:Any})

Removes tags from a detector, dataset, or alert.

# Arguments
- `resource_arn`: The resource's Amazon Resource Name (ARN).
- `tag_keys`: Keys to remove from the resource's tags.

"""
untag_resource(resourceArn, tagKeys; aws_config::AbstractAWSConfig=global_aws_config()) = lookoutmetrics("DELETE", "/tags/$(resourceArn)", Dict{String, Any}("tagKeys"=>tagKeys); aws_config=aws_config)
untag_resource(resourceArn, tagKeys, params::AbstractDict{String}; aws_config::AbstractAWSConfig=global_aws_config()) = lookoutmetrics("DELETE", "/tags/$(resourceArn)", Dict{String, Any}(mergewith(_merge, Dict{String, Any}("tagKeys"=>tagKeys), params)); aws_config=aws_config)

"""
    update_anomaly_detector(anomaly_detector_arn)
    update_anomaly_detector(anomaly_detector_arn, params::Dict{String,<:Any})

Updates a detector. After activation, you can only change a detector's ingestion delay and
description.

# Arguments
- `anomaly_detector_arn`: The ARN of the detector to update.

# Optional Parameters
Optional parameters can be passed as a `params::Dict{String,<:Any}`. Valid keys are:
- `"AnomalyDetectorConfig"`: Contains information about the configuration to which the
  detector will be updated.
- `"AnomalyDetectorDescription"`: The updated detector description.
- `"KmsKeyArn"`: The Amazon Resource Name (ARN) of an AWS KMS encryption key.
"""
update_anomaly_detector(AnomalyDetectorArn; aws_config::AbstractAWSConfig=global_aws_config()) = lookoutmetrics("POST", "/UpdateAnomalyDetector", Dict{String, Any}("AnomalyDetectorArn"=>AnomalyDetectorArn); aws_config=aws_config)
update_anomaly_detector(AnomalyDetectorArn, params::AbstractDict{String}; aws_config::AbstractAWSConfig=global_aws_config()) = lookoutmetrics("POST", "/UpdateAnomalyDetector", Dict{String, Any}(mergewith(_merge, Dict{String, Any}("AnomalyDetectorArn"=>AnomalyDetectorArn), params)); aws_config=aws_config)

"""
    update_metric_set(metric_set_arn)
    update_metric_set(metric_set_arn, params::Dict{String,<:Any})

Updates a dataset.

# Arguments
- `metric_set_arn`: The ARN of the dataset to update.

# Optional Parameters
Optional parameters can be passed as a `params::Dict{String,<:Any}`. Valid keys are:
- `"DimensionList"`: The dimension list.
- `"MetricList"`: The metric list.
- `"MetricSetDescription"`: The dataset's description.
- `"MetricSetFrequency"`: The dataset's interval.
- `"MetricSource"`:
- `"Offset"`: After an interval ends, the amount of time that the detector waits before
  importing data.
- `"TimestampColumn"`: The timestamp column.
"""
update_metric_set(MetricSetArn; aws_config::AbstractAWSConfig=global_aws_config()) = lookoutmetrics("POST", "/UpdateMetricSet", Dict{String, Any}("MetricSetArn"=>MetricSetArn); aws_config=aws_config)
update_metric_set(MetricSetArn, params::AbstractDict{String}; aws_config::AbstractAWSConfig=global_aws_config()) = lookoutmetrics("POST", "/UpdateMetricSet", Dict{String, Any}(mergewith(_merge, Dict{String, Any}("MetricSetArn"=>MetricSetArn), params)); aws_config=aws_config)


text:
module EnsembleKalmanProcessModule

using ..ParameterDistributionStorage
using ..DataStorage

using Random
using Statistics
using Distributions
using LinearAlgebra
using DocStringExtensions

export EnsembleKalmanProcess
export get_u, get_g
export get_u_prior, get_u_final, get_u_mean_final, get_g_final, get_N_iterations, get_error
export compute_error!
export update_ensemble!

abstract type Process end
#specific Processes and their exports are included after the general definitions

## begin general constructor and function definitions

"""
    EnsembleKalmanProcess{FT<:AbstractFloat, IT<:Int}

Structure that is used in Ensemble Kalman processes

#Fields
$(DocStringExtensions.FIELDS)
"""
struct EnsembleKalmanProcess{FT<:AbstractFloat, IT<:Int, P<:Process}
    "Array of stores for parameters (u), each of size [N_par × N_ens]"
    u::Array{DataContainer{FT}}
    "vector of the observed vector size [N_obs]"
    obs_mean::Vector{FT}
    "covariance matrix of the observational noise, of size [N_obs × N_obs]"
    obs_noise_cov::Array{FT, 2}
    "ensemble size"
    N_ens::IT
    "Array of stores for forward model outputs, each of size  [N_obs × N_ens]"
    g::Array{DataContainer{FT}}
    "vector of errors"
    err::Vector{FT}
    "vector of timesteps used in each EK iteration"
    Δt::Vector{FT}
    "the particular EK process (`Inversion` or `Sampler` or `Unscented`)"
    process::P
end

# outer constructors
function EnsembleKalmanProcess(params::Array{FT, 2},
                               obs_mean,
                               obs_noise_cov::Array{FT, 2},
                               process::P;
                               Δt=FT(1)) where {FT<:AbstractFloat, P<:Process}

    #initial parameters stored as columns
    init_params=DataContainer(params, data_are_columns=true)
    # ensemble size
    N_ens = size(init_params,2) #stored with data as columns
    IT = typeof(N_ens)
    #store for model evaluations
    g=[]
    # error store
    err = FT[]
    # timestep store
    Δt = Array([Δt])

    EnsembleKalmanProcess{FT, IT, P}([init_params], obs_mean, obs_noise_cov, N_ens, g,
                                     err, Δt, process)
end

"""
    get_u(ekp::EnsembleKalmanProcess, iteration::IT; return_array=true) where {IT <: Integer}

Get for the EKI iteration. Returns a DataContainer object unless array is specified.
"""
function get_u(ekp::EnsembleKalmanProcess, iteration::IT; return_array=true) where {IT <: Integer}
    return  return_array ? get_data(ekp.u[iteration]) : ekp.u[iteration]
end

"""
    get_g(ekp::EnsembleKalmanProcess, iteration::IT; return_array=true) where {IT <: Integer}

Get for the EKI iteration. Returns a DataContainer object unless array is specified.
"""
function get_g(ekp::EnsembleKalmanProcess, iteration::IT; return_array=true) where {IT <: Integer}
    return return_array ? get_data(ekp.g[iteration]) : ekp.g[iteration]
end

"""
    get_u(ekp::EnsembleKalmanProcess; return_array=true)

Get for the EKI iteration. Returns a DataContainer object unless array is specified.
"""
function get_u(ekp::EnsembleKalmanProcess; return_array=true) where {IT <: Integer}
    N_stored_u = get_N_iterations(ekp)+1
    return [get_u(ekp, it, return_array=return_array) for it in 1:N_stored_u]
end

"""
    get_g(ekp::EnsembleKalmanProcess; return_array=true)

Get for the EKI iteration. Returns a DataContainer object unless array is specified.
"""
function get_g(ekp::EnsembleKalmanProcess; return_array=true) where {IT <: Integer}
    N_stored_g = get_N_iterations(ekp)
    return [get_g(ekp, it, return_array=return_array) for it in 1:N_stored_g]
end

"""
    get_u_final(ekp::EnsembleKalmanProcess, return_array=true)

Get the final or prior iteration of parameters or model ouputs, returns a DataContainer Object if return_array is false.
"""
function get_u_final(ekp::EnsembleKalmanProcess; return_array=true)
    return return_array ? get_u(ekp,size(ekp.u)[1]) : ekp.u[end]
end

"""
    get_u_prior(ekp::EnsembleKalmanProcess, return_array=true)

Get the final or prior iteration of parameters or model ouputs, returns a DataContainer Object if return_array is false.
"""

function get_u_prior(ekp::EnsembleKalmanProcess; return_array=true)
    return return_array ? get_u(ekp,1) : ekp.u[1]
end

"""
    get_g_final(ekp::EnsembleKalmanProcess, return_array=true)

Get the final or prior iteration of parameters or model ouputs, returns a DataContainer Object if return_array is false.
"""

function get_g_final(ekp::EnsembleKalmanProcess; return_array=true)
    return return_array ? get_g(ekp,size(ekp.g)[1]) : ekp.g[end]
end

"""
    get_N_iterations(ekp::EnsembleKalmanProcess

get number of times update has been called (equals size(g), or size(u)-1)
"""
function get_N_iterations(ekp::EnsembleKalmanProcess)
    return size(ekp.u)[1] - 1
end

"""
    construct_initial_ensemble(prior::ParameterDistribution, N_ens::IT; rng_seed=42) where {IT<:Int}

Construct the initial parameters, by sampling N_ens samples from specified
prior distribution. Returned with parameters as columns
"""
function construct_initial_ensemble(prior::ParameterDistribution, N_ens::IT; rng_seed=42) where {IT<:Int}
    # Ensuring reproducibility of the sampled parameter values
    Random.seed!(rng_seed)
    parameters = sample_distribution(prior, N_ens) #of size [dim(param space) N_ens]
    return parameters
end

function compute_error!(ekp::EnsembleKalmanProcess)
    mean_g = dropdims(mean(get_g_final(ekp), dims=2), dims=2)
    diff = ekp.obs_mean - mean_g
    X = ekp.obs_noise_cov \ diff # diff: column vector
    newerr = dot(diff, X)
    push!(ekp.err, newerr)
end

function get_error(ekp::EnsembleKalmanProcess)
    return ekp.err
end

## include the different types of Processes and their exports:

# struct Inversion
export Inversion
export find_ekp_stepsize
include("EnsembleKalmanInversion.jl")

# struct Sampler
export Sampler
include("EnsembleKalmanSampler.jl")

# struct Unscented
export Unscented
export construct_initial_ensemble
export Gaussian_2d
include("UnscentedKalmanInversion.jl")

end # module


text:
# by lazarusA # HIDE
#using GLMakie, Random # HIDE
using CairoMakie, Random
CairoMakie.activate!() #HIDE
let
    Random.seed!(13)
    n = 3000
    data = randn(n)
    fig = Figure(resolution = (700, 450), font = "sans", fontsize = 20)
    ax1 = Axis(fig, xlabel = "value" )
    hist!(ax1,data,normalization=:pdf,color = (:red,0.5), label = "hist & pdf")
    density!(ax1, data, color =(:grey, 0.5), label = "density!", strokewidth=1)
    axislegend(ax1, position = :rt)
    fig[1,1] = ax1
    fig
    save(joinpath(@__DIR__, "output", "histDensity.png"), current_figure(), px_per_unit = 2.0) # HIDE
end
using Pkg # HIDE
Pkg.status(["CairoMakie", "Random"]) # HIDE

max_stars_repo_name: MisterBiggs/BeautifulMakie

max_issues_repo_name: MisterBiggs/BeautifulMakie

max_forks_repo_name: MisterBiggs/BeautifulMakie


text:
using NDPrototype
using LightGraphs
using Random
using OrdinaryDiffEq

#############
include("../benchmark/benchmarks.jl")
include("../benchmark/benchmark_utils.jl")

# diffusion test
begin
    N = 100000
    g = watts_strogatz(N, 3, 0.8, seed=1)
    edge = diffusion_dedge()
    vertex = diffusion_vertex()
    nd1 = Network(g, vertex, edge, execution=:seq);
    nd2 = Network(g, vertex, edge, execution=:threaded);
    x0 = randn(dim(nd1));
    dx1 = zero(x0);
    dx2 = zero(x0);
    nd1(dx1, x0, nothing, 0.0)
    nd2(dx2, x0, nothing, 0.0)
    @test dx1 ≈ dx2

    @btime $nd1($dx1, $x0, nothing, 0.0)
    @btime $nd2($dx2, $x0, nothing, 0.0)

    prob1 = ODEProblem(nd1, x0, (0.0, 1.0), nothing)
    prob2 = ODEProblem(nd2, x0, (0.0, 1.0), nothing)
    @time sol1a = solve(prob1, Tsit5());
    @time sol1b = solve(prob1, Rodas4());
    @time sol1c = solve(prob1, Rosenbrock23());
    @time sol2a = solve(prob2, Tsit5());
    @time sol2b = solve(prob2, Rodas4());
end

begin
    N = 100_000
    (p, v, e, g) = heterogeneous(N)
    nd1 = Network(g, v, e, execution=:seq)
    (p, v, e, g) = heterogeneous(N)
    nd2 = Network(g, v, e, execution=:threaded)
    x0 = randn(dim(nd1))
    dx1 = zeros(dim(nd1))
    dx2 = zeros(dim(nd1))
    @time nd1(dx1, x0, p, 0.0) # call to init caches, we don't want to benchmark this
    @time nd1(dx2, x0, p, 0.0) # call to init caches, we don't want to benchmark this
    @test dx1 ≈ dx2

    @time nd2(dx1, x0, p, 0.0) # call to init caches, we don't want to benchmark this
    @time nd2(dx2, x0, p, 0.0) # call to init caches, we don't want to benchmark this
    @test dx1 ≈ dx2

    @btime $nd1($dx1, $x0, $p, 0.0)
  # 509.375 μs (6 allocations: 320 bytes)
    @btime $nd2($dx2, $x0, $p, 0.0)
  # 211.208 μs (222 allocations: 17.84 KiB)

    N = 1000
    g = watts_strogatz(N, Int(N/2), 0.0, seed=1)
    gc = NDPrototype.ColoredGraph(g)
    for e in LightGraphs.edges(gc)
        c = NDPrototype.pickfree(gc, e.src, e.dst)
        NDPrototype.setcolor!(gc, e, c)
        println(e)
    end
    gc
end