text:
# using ExtendableSparse
# using ExtendableGrids

function get_testgrid(::Type{<:Edge1D})
    X=collect(0:0.05:1)
    simplexgrid(X)
end

function get_testgrid(::Type{<:Triangle2D})
    X=collect(0:0.05:1)
    Y=collect(0:0.05:1)
    simplexgrid(X,Y)
end

function get_testgrid(::Type{<:Tetrahedron3D})
    X=collect(0:0.1:1)
    Y=collect(0:0.1:1)
    Z=collect(0:0.1:1)
    simplexgrid(X,Y,Z)
end

function get_testgrid(::Type{Parallelogram2D})
    return uniform_refine(grid_unitsquare(Parallelogram2D),2)
end

function get_testgrid(::Type{Parallelepiped3D})
    return uniform_refine(grid_unitcube(Parallelepiped3D),2)
end

function get_testgrid(::Type{Triangle2D},::Type{Parallelogram2D})
    return uniform_refine(grid_unitsquare_mixedgeometries(),1)
end

function check_enumeration_consistency(G::Type{<:AbstractElementGeometry}, CellItems, ItemNodes, item_enum_rule)

    @info "Checking enumeration consistency $CellItems <> $ItemNodes with $item_enum_rule for geometry = $G..."

    xgrid = get_testgrid(G)

    xCellNodes = xgrid[CellNodes]
    xCellItems = xgrid[CellItems]
    xItemNodes = xgrid[ItemNodes]
    item_rule = item_enum_rule(G)
    nitems_per_cell = max_num_targets_per_source(xCellItems)
    nnodes_per_item = max_num_targets_per_source(xItemNodes)
    ncells = num_sources(xCellNodes)

    @assert nitems_per_cell == size(item_rule,2)
    @assert nnodes_per_item == size(item_rule,1)

    item::Int = 0
    everything_is_consistent = true
    itemnodes_found = zeros(Bool, nnodes_per_item)
    for cell = 1 : ncells
        for i = 1 : nitems_per_cell
            item = xCellItems[i,cell]
            fill!(itemnodes_found,false)
            for n = 1 : nnodes_per_item
                for k = 1 : nnodes_per_item
                    if xItemNodes[n,item] == xCellNodes[item_rule[k,i],cell]
                        itemnodes_found[n] = true
                        break
                    end
                end
            end
            if prod(itemnodes_found) != true
                everything_is_consistent = false
            end
        end
    end

    return everything_is_consistent
end

function check_cellfinder(xgrid)
    EG = xgrid[UniqueCellGeometries]
    @info "Testing CellFinder for geometries=$EG..."
    xCoordinates = xgrid[Coordinates]
    xCellNodes = xgrid[CellNodes]
    edim = dim_element(EG[1])
    CF::CellFinder{Float64,Int32} = CellFinder(xgrid)
    cell_to_find = num_sources(xCellNodes) - 1

    # compute midpoint of cell_to_find
    x_source = zeros(Float64,edim)
    for j = 1 : edim, k = 1 : num_targets(xCellNodes,cell_to_find)
        x_source[j] += xCoordinates[j,xCellNodes[k,cell_to_find]] + 1e-6
    end
    x_source ./= num_targets(xCellNodes,cell_to_find)

    # find cell by local strategy
    xref = zeros(Float64,edim+1)
    cell = gFindLocal!(xref, CF, x_source; icellstart = 1)

    # check xref
    x = zeros(Float64,edim)
    L2G = L2GTransformer(xgrid[CellGeometries][cell], xgrid, ON_CELLS)
    update_trafo!(L2G,cell)
    eval_trafo!(x,L2G,xref)

    @info "... found x=$x in cell = $cell by local search (and had to find x=$x_source in cell=$cell_to_find)"

    @assert cell == cell_to_find

    # find cell again by brute force
    cell = gFindBruteForce!(xref, CF, x_source)

    # check xref
    x = zeros(Float64,edim)
    update_trafo!(L2G,cell)
    eval_trafo!(x,L2G,xref)

    @info "... found x=$x in cell = $cell by brute force (and had to find  x=$x_source in cell=$cell_to_find)"

    @assert cell == cell_to_find

    return (cell == cell_to_find) && (sqrt(sum((x - x_source).^2)) < 1e-14)
end

function run_grid_tests()
    # test FACE enumerations
    @test check_enumeration_consistency(Edge1D, CellFaces, FaceNodes, local_cellfacenodes)
    @test check_enumeration_consistency(Triangle2D, CellFaces, FaceNodes, local_cellfacenodes)
    @test check_enumeration_consistency(Parallelogram2D, CellFaces, FaceNodes, local_cellfacenodes)
    @test check_enumeration_consistency(Tetrahedron3D, CellFaces, FaceNodes, local_cellfacenodes)
    @test check_enumeration_consistency(Parallelepiped3D, CellFaces, FaceNodes, local_cellfacenodes)

    # test EDGE enumerations
    @test check_enumeration_consistency(Tetrahedron3D, CellEdges, EdgeNodes, local_celledgenodes)
    @test check_enumeration_consistency(Parallelepiped3D, CellEdges, EdgeNodes, local_celledgenodes)

    # todo: check FaceEdges, CellFaceSigns, CellFaceOrientations
    # todo: FaceCells, EdgeCells

    # test CellFinder
    @test check_cellfinder(get_testgrid(Edge1D))
    @test check_cellfinder(get_testgrid(Triangle2D))
    @test check_cellfinder(get_testgrid(Parallelogram2D))
    @test check_cellfinder(get_testgrid(Triangle2D,Parallelogram2D))
    @test check_cellfinder(get_testgrid(Tetrahedron3D))
    @test check_cellfinder(get_testgrid(Parallelepiped3D))
end


text:
module TestConstraintToInterval
# These tests are mostly copies of the flip_sign.jl tests for GreaterToLess

using Test

using MathOptInterface
const MOI = MathOptInterface

function runtests()
    for name in names(@__MODULE__; all = true)
        if startswith("$(name)", "test_")
            @testset "$(name)" begin
                getfield(@__MODULE__, name)()
            end
        end
    end
    return
end

include("../utilities.jl")

function test_GreaterToInterval()
    mock = MOI.Utilities.MockOptimizer(
        MOI.Utilities.UniversalFallback(MOI.Utilities.Model{Float64}()),
    )
    config = MOI.Test.Config()
    bridged_mock = MOI.Bridges.Constraint.GreaterToInterval{Float64}(mock)
    MOI.Test.runtests(
        bridged_mock,
        config,
        include = [
            "test_basic_$(F)_GreaterThanThan" for F in
            ["VariableIndex", "ScalarAffineFunction", "ScalarQuadraticFunction"]
        ],
    )
    MOI.empty!(bridged_mock)
    MOI.Utilities.set_mock_optimize!(
        mock,
        (mock::MOI.Utilities.MockOptimizer) ->
            MOI.Utilities.mock_optimize!(mock, [0.0, 0.0]),
        (mock::MOI.Utilities.MockOptimizer) ->
            MOI.Utilities.mock_optimize!(mock, [100.0, 0.0]),
        (mock::MOI.Utilities.MockOptimizer) ->
            MOI.Utilities.mock_optimize!(mock, [100.0, -100.0]),
    )
    MOI.Test.test_linear_modify_GreaterThan_and_LessThan_constraints(
        bridged_mock,
        config,
    )
    ci = first(
        MOI.get(
            bridged_mock,
            MOI.ListOfConstraintIndices{
                MOI.ScalarAffineFunction{Float64},
                MOI.GreaterThan{Float64},
            }(),
        ),
    )
    for attr in [MOI.ConstraintPrimalStart(), MOI.ConstraintDualStart()]
        @test MOI.supports(bridged_mock, attr, typeof(ci))
        MOI.set(bridged_mock, attr, ci, 2.0)
        @test MOI.get(bridged_mock, attr, ci) ≈ 2.0
    end
    _test_delete_bridge(
        bridged_mock,
        ci,
        2,
        ((MOI.ScalarAffineFunction{Float64}, MOI.Interval{Float64}, 0),),
    )
    return
end

function test_LessToInterval()
    mock = MOI.Utilities.MockOptimizer(
        MOI.Utilities.UniversalFallback(MOI.Utilities.Model{Float64}()),
    )
    config = MOI.Test.Config()
    bridged_mock = MOI.Bridges.Constraint.LessToInterval{Float64}(mock)
    MOI.Test.runtests(
        bridged_mock,
        config,
        include = [
            "test_basic_$(F)_LessThan" for F in
            ["VariableIndex", "ScalarAffineFunction", "ScalarQuadraticFunction"]
        ],
    )
    MOI.empty!(bridged_mock)
    MOI.Utilities.set_mock_optimize!(
        mock,
        (mock::MOI.Utilities.MockOptimizer) -> MOI.Utilities.mock_optimize!(
            mock,
            MOI.OPTIMAL,
            (MOI.FEASIBLE_POINT, [1.0]),
            MOI.FEASIBLE_POINT,
            (MOI.ScalarAffineFunction{Float64}, MOI.Interval{Float64}) =>
                [-1.0],
        ),
        (mock::MOI.Utilities.MockOptimizer) -> MOI.Utilities.mock_optimize!(
            mock,
            MOI.OPTIMAL,
            (MOI.FEASIBLE_POINT, [2.0]),
            MOI.FEASIBLE_POINT,
            (MOI.ScalarAffineFunction{Float64}, MOI.Interval{Float64}) =>
                [-1.0],
        ),
    )
    MOI.Test.test_modification_set_scalaraffine_lessthan(bridged_mock, config)
    MOI.empty!(bridged_mock)
    MOI.Utilities.set_mock_optimize!(
        mock,
        (mock::MOI.Utilities.MockOptimizer) -> MOI.Utilities.mock_optimize!(
            mock,
            MOI.OPTIMAL,
            (MOI.FEASIBLE_POINT, [1.0]),
            MOI.FEASIBLE_POINT,
            (MOI.ScalarAffineFunction{Float64}, MOI.Interval{Float64}) =>
                [-1.0],
        ),
        (mock::MOI.Utilities.MockOptimizer) -> MOI.Utilities.mock_optimize!(
            mock,
            MOI.OPTIMAL,
            (MOI.FEASIBLE_POINT, [0.5]),
            MOI.FEASIBLE_POINT,
            (MOI.ScalarAffineFunction{Float64}, MOI.Interval{Float64}) =>
                [-0.5],
        ),
    )
    MOI.Test.test_modification_coef_scalaraffine_lessthan(bridged_mock, config)
    ci = first(
        MOI.get(
            bridged_mock,
            MOI.ListOfConstraintIndices{
                MOI.ScalarAffineFunction{Float64},
                MOI.LessThan{Float64},
            }(),
        ),
    )
    for attr in [MOI.ConstraintPrimalStart(), MOI.ConstraintDualStart()]
        @test MOI.supports(bridged_mock, attr, typeof(ci))
        MOI.set(bridged_mock, attr, ci, 2.0)
        @test MOI.get(bridged_mock, attr, ci) ≈ 2.0
    end
    _test_delete_bridge(
        bridged_mock,
        ci,
        1,
        ((MOI.ScalarAffineFunction{Float64}, MOI.Interval{Float64}, 0),),
    )
    return
end

# Define a dummy optimizer that only supports intervals
# and use it in the below unmocked test
mutable struct Optimizer <: MOI.AbstractOptimizer end

MOI.get(model::Optimizer, ::MOI.SolverName) = "OnlyIntervalOptimizer"
MOI.get(model::Optimizer, ::MOI.SolverVersion) = "v0.1.0"

function MOI.supports_constraint(
    ::Optimizer,
    ::Type{MOI.ScalarAffineFunction{Float64}},
    ::Type{MOI.Interval{Float64}},
)
    return true
end

function test_GreaterOrLessToInterval_unmocked()
    # model supports Interval but not LessThan or GreaterThan
    model = Optimizer()
    @test MOI.supports_constraint(
        model,
        MOI.ScalarAffineFunction{Float64},
        MOI.Interval{Float64},
    )
    @test !MOI.supports_constraint(
        model,
        MOI.ScalarAffineFunction{Float64},
        MOI.LessThan{Float64},
    )
    @test !MOI.supports_constraint(
        model,
        MOI.ScalarAffineFunction{Float64},
        MOI.GreaterThan{Float64},
    )

    # bridged model supports all
    bridged = MOI.Bridges.Constraint.GreaterToInterval{Float64}(
        MOI.Bridges.Constraint.LessToInterval{Float64}(model),
    )
    @test MOI.supports_constraint(
        bridged,
        MOI.ScalarAffineFunction{Float64},
        MOI.Interval{Float64},
    )
    @test MOI.supports_constraint(
        bridged,
        MOI.ScalarAffineFunction{Float64},
        MOI.LessThan{Float64},
    )
    @test MOI.supports_constraint(
        bridged,
        MOI.ScalarAffineFunction{Float64},
        MOI.GreaterThan{Float64},
    )

    # bridged model with Bridges.full_bridge_optimizer
    bridged2 = MOI.Bridges.full_bridge_optimizer(model, Float64)
    @test MOI.supports_constraint(
        bridged2,
        MOI.ScalarAffineFunction{Float64},
        MOI.Interval{Float64},
    )
    @test MOI.supports_constraint(
        bridged2,
        MOI.ScalarAffineFunction{Float64},
        MOI.LessThan{Float64},
    )
    @test MOI.supports_constraint(
        bridged2,
        MOI.ScalarAffineFunction{Float64},
        MOI.GreaterThan{Float64},
    )
    return
end

end  # module

TestConstraintToInterval.runtests()

ext: jl


text:
#!/usr/bin/env julia

using CairoMakie

CairoMakie.activate!(type = "png")

x = range(0, 10, length=100)
y = sin.(x)
scene = lines(x, y)
save("plot.png",scene)

max_stars_repo_name: jdurbin/sandbox

max_issues_repo_name: jdurbin/sandbox

max_forks_repo_name: jdurbin/sandbox


text:
export nframes,
    frame_type,
    select_last_dim,
    select_last_frame,
    consecutive_view,
    find_all_max,
    huber_loss,
    huber_loss_unreduced,
    discount_rewards,
    discount_rewards!,
    discount_rewards_reduced,
    generalized_advantage_estimation,
    generalized_advantage_estimation!,
    logitcrossentropy_unreduced,
    flatten_batch,
    unflatten_batch

using StatsBase

nframes(a::AbstractArray{T,N}) where {T,N} = size(a, N)
frame_type(::Array{T,N}) where {T,N} = Array{T,N - 1}
frame_type(::Vector{T}) where {T} = T

select_last_dim(xs::AbstractArray{T,N}, inds) where {T,N} =
    @views xs[ntuple(_ -> (:), N - 1)..., inds]

select_last_frame(xs::AbstractArray{T,N}) where {T,N} = select_last_dim(xs, size(xs, N))

"""
    flatten_batch(x::AbstractArray)

Merge the last two dimension.

# Example

```julia-repl
julia> x = reshape(1:12, 2, 2, 3)
2×2×3 reshape(::UnitRange{Int64}, 2, 2, 3) with eltype Int64:
[:, :, 1] =
 1  3
 2  4

[:, :, 2] =
 5  7
 6  8

[:, :, 3] =
  9  11
 10  12

julia> flatten_batch(x)
2×6 reshape(::UnitRange{Int64}, 2, 6) with eltype Int64:
 1  3  5  7   9  11
 2  4  6  8  10  12
```
"""
flatten_batch(x::AbstractArray) =
    reshape(x, (size(x) |> reverse |> Base.tail |> Base.tail |> reverse)..., :)  # much faster than  `reshape(x, size(x)[1:end-2]..., :)`

unflatten_batch(x::AbstractArray, i::Int...) =
    reshape(x, (size(x) |> reverse |> Base.tail |> reverse)..., i...)

consecutive_view(
    cb::AbstractArray,
    inds::Vector{Int};
    n_stack = nothing,
    n_horizon = nothing,
) = consecutive_view(cb, inds, n_stack, n_horizon)
consecutive_view(cb::AbstractArray, inds::Vector{Int}, ::Nothing, ::Nothing) =
    select_last_dim(cb, inds)
consecutive_view(cb::AbstractArray, inds::Vector{Int}, n_stack::Int, ::Nothing) =
    select_last_dim(
        cb,
        reshape([i for x in inds for i in x-n_stack+1:x], n_stack, length(inds)),
    )
consecutive_view(cb::AbstractArray, inds::Vector{Int}, ::Nothing, n_horizon::Int) =
    select_last_dim(
        cb,
        reshape([i for x in inds for i in x:x+n_horizon-1], n_horizon, length(inds)),
    )
consecutive_view(cb::AbstractArray, inds::Vector{Int}, n_stack::Int, n_horizon::Int) =
    select_last_dim(
        cb,
        reshape(
            [j for x in inds for i in x:x+n_horizon-1 for j in i-n_stack+1:i],
            n_stack,
            n_horizon,
            length(inds),
        ),
    )

function find_all_max(x)
    v = maximum(x)
    v, findall(==(v), x)
end

function find_all_max(x, mask::AbstractVector{Bool})
    v = maximum(view(x, mask))
    v, [k for (m, k) in zip(mask, keys(x)) if m && x[k] == v]
end

# !!! watch https://github.com/JuliaLang/julia/pull/35316#issuecomment-622629895
Base.findmax(f, domain) = mapfoldl(x -> (f(x), x), _rf_findmax, domain)
_rf_findmax((fm, m), (fx, x)) = isless(fm, fx) ? (fx, x) : (fm, m)

# !!! type piracy
Base.findmax(A::AbstractVector, mask::AbstractVector{Bool}) =
    findmax(i -> A[i], view(keys(A), mask))

function logitcrossentropy_unreduced(logŷ::AbstractVecOrMat, y::AbstractVecOrMat)
    return vec(-sum(y .* logsoftmax(logŷ), dims = 1))
end

"""
    huber_loss_unreduced(labels, predictions; δ = 1.0f0)

Similar to [`huber_loss`](@ref), but it doesn't do the `mean` operation in the last step.
"""
function huber_loss_unreduced(labels, predictions; δ = 1.0f0)
    abs_error = abs.(predictions .- labels)
    quadratic = min.(abs_error, δ)
    linear = abs_error .- quadratic
    0.5f0 .* quadratic .* quadratic .+ δ .* linear
end

"""
    huber_loss(labels, predictions; δ = 1.0f0)

See [Huber loss](https://en.wikipedia.org/wiki/Huber_loss)
"""
huber_loss(labels, predictions; δ = 1.0f0) =
    huber_loss_unreduced(labels, predictions; δ = δ) |> mean

const VectorOrMatrix = Union{AbstractMatrix,AbstractVector}

"""
    discount_rewards(rewards::VectorOrMatrix, γ::Number;kwargs...)

Calculate the gain started from the current step with discount rate of `γ`.
`rewards` can be a matrix.

# Keyword argments

- `dims=:`, if `rewards` is a `Matrix`, then `dims` can only be `1` or `2`.
- `terminal=nothing`, specify if each reward follows by a terminal. `nothing` means the game is not terminated yet. If `terminal` is provided, then the size must be the same with `rewards`.
- `init=nothing`, `init` can be used to provide the the reward estimation of the last state.

# Example
"""
function discount_rewards(rewards::VectorOrMatrix, γ::T; kwargs...) where {T<:Number}
    res = similar(rewards, promote_type(eltype(rewards), T))
    discount_rewards!(res, rewards, γ; kwargs...)
    res
end

discount_rewards!(new_rewards, rewards, γ; terminal = nothing, init = nothing, dims = :) =
    _discount_rewards!(new_rewards, rewards, γ, terminal, init, dims)

function _discount_rewards!(
    new_rewards::AbstractMatrix,
    rewards::AbstractMatrix,
    γ,
    terminal::Nothing,
    init::Nothing,
    dims::Int,
)
    dims = ndims(rewards) - dims + 1
    for (r′, r) in zip(eachslice(new_rewards, dims = dims), eachslice(rewards, dims = dims))
        _discount_rewards!(r′, r, γ, nothing, nothing)
    end
end

function _discount_rewards!(
    new_rewards::AbstractMatrix,
    rewards::AbstractMatrix,
    γ,
    terminal::Nothing,
    init,
    dims::Int,
)
    dims = ndims(rewards) - dims + 1
    for (i, (r′, r)) in
        enumerate(zip(eachslice(new_rewards, dims = dims), eachslice(rewards, dims = dims)))
        _discount_rewards!(r′, r, γ, nothing, init[i])
    end
end

function _discount_rewards!(
    new_rewards::AbstractMatrix,
    rewards::AbstractMatrix,
    γ,
    terminal,
    init::Nothing,
    dims::Int,
)
    dims = ndims(rewards) - dims + 1
    for (r′, r, t) in zip(
        eachslice(new_rewards, dims = dims),
        eachslice(rewards, dims = dims),
        eachslice(terminal, dims = dims),
    )
        _discount_rewards!(r′, r, γ, t, nothing)
    end
end

function _discount_rewards!(
    new_rewards::AbstractMatrix,
    rewards::AbstractMatrix,
    γ,
    terminal,
    init,
    dims::Int,
)
    dims = ndims(rewards) - dims + 1
    for (i, (r′, r, t)) in enumerate(zip(
        eachslice(new_rewards, dims = dims),
        eachslice(rewards, dims = dims),
        eachslice(terminal, dims = dims),
    ))
        _discount_rewards!(r′, r, γ, t, init[i])
    end
end

_discount_rewards!(
    new_rewards::AbstractVector,
    rewards::AbstractVector,
    γ,
    terminal,
    init,
    dims::Colon,
) = _discount_rewards!(new_rewards, rewards, γ, terminal, init)

"assuming rewards and new_rewards are Vector"
_discount_rewards!(new_rewards, rewards, γ, terminal, init::Nothing) =
    _discount_rewards!(new_rewards, rewards, γ, terminal, zero(eltype(new_rewards)))

function _discount_rewards!(new_rewards, rewards, γ, terminal, init)
    gain = init
    for i in length(rewards):-1:1
        is_continue = isnothing(terminal) ? true : (!terminal[i])
        gain = rewards[i] + γ * gain * is_continue
        new_rewards[i] = gain
    end
    new_rewards
end

discount_rewards_reduced(rewards::AbstractVector, γ; terminal = nothing, init = nothing) =
    _discount_rewards_reduced(rewards, γ, terminal, init)

function discount_rewards_reduced(
    rewards::AbstractMatrix,
    γ::T;
    terminal = nothing,
    init = nothing,
    dims,
) where {T<:Number}
    dims = ndims(rewards) - dims + 1
    res = Array{promote_type(eltype(rewards), T)}(undef, size(rewards, dims))
    _discount_rewards_reduced!(res, rewards, γ, terminal, init, dims)
    res
end

_discount_rewards_reduced(rewards, γ, terminal, init::Nothing) =
    _discount_rewards_reduced(rewards, γ, terminal, zero(eltype(rewards)))

function _discount_rewards_reduced(rewards, γ, terminal, init)
    gain = init
    for i in length(rewards):-1:1
        is_continue = isnothing(terminal) ? true : (!terminal[i])
        gain = rewards[i] + γ * gain * is_continue
    end
    gain
end

discount_rewards_reduced!(
    reduced_rewards::AbstractVector,
    rewards::AbstractMatrix,
    γ;
    terminal = nothing,
    init = nothing,
    dims,
) = _discount_rewards_reduced!(reduced_rewards, rewards, γ, terminal, init, dims)

function _discount_rewards_reduced!(
    reduced_rewards,
    rewards,
    γ,
    terminal::Nothing,
    init::Nothing,
    dims::Int,
)
    for (i, r) in enumerate(eachslice(rewards, dims = dims))
        reduced_rewards[i] = _discount_rewards_reduced(r, γ, nothing, nothing)
    end
end

function _discount_rewards_reduced!(
    reduced_rewards,
    rewards,
    γ,
    terminal::Nothing,
    init,
    dims::Int,
)
    for (i, r) in enumerate(eachslice(rewards, dims = dims))
        reduced_rewards[i] = _discount_rewards_reduced(r, γ, nothing, init[i])
    end
end

function _discount_rewards_reduced!(
    reduced_rewards,
    rewards,
    γ,
    terminal,
    init::Nothing,
    dims::Int,
)
    for (i, (r, t)) in
        enumerate(zip(eachslice(rewards, dims = dims), eachslice(terminal, dims = dims)))
        reduced_rewards[i] = _discount_rewards_reduced(r, γ, t, nothing)
    end
end

function _discount_rewards_reduced!(reduced_rewards, rewards, γ, terminal, init, dims::Int)
    for (i, (r, t)) in
        enumerate(zip(eachslice(rewards, dims = dims), eachslice(terminal, dims = dims)))
        reduced_rewards[i] = _discount_rewards_reduced(r, γ, t, init[i])
    end
end

"""
    generalized_advantage_estimation(rewards::VectorOrMatrix, values::VectorOrMatrix, γ::Number, λ::Number;kwargs...)

Calculate the generalized advantage estimate started from the current step with discount rate of `γ` and a lambda for GAE-Lambda of 'λ'.
`rewards` and 'values' can be a matrix.

# Keyword argments

- `dims=:`, if `rewards` is a `Matrix`, then `dims` can only be `1` or `2`.
- `terminal=nothing`, specify if each reward follows by a terminal. `nothing` means the game is not terminated yet. If `terminal` is provided, then the size must be the same with `rewards`.

# Example
"""
function generalized_advantage_estimation(
    rewards::VectorOrMatrix,
    values::VectorOrMatrix,
    γ::T,
    λ::T;
    kwargs...,
) where {T<:Number}
    res = similar(rewards, promote_type(eltype(rewards), T))
    generalized_advantage_estimation!(res, rewards, values, γ, λ; kwargs...)
    res
end

generalized_advantage_estimation!(
    advantages,
    rewards,
    values,
    γ,
    λ;
    terminal = nothing,
    dims = :,
) = _generalized_advantage_estimation!(advantages, rewards, values, γ, λ, terminal, dims)

function _generalized_advantage_estimation!(
    advantages::AbstractMatrix,
    rewards::AbstractMatrix,
    values::AbstractMatrix,
    γ,
    λ,
    terminal::Nothing,
    dims::Int,
)
    dims = ndims(rewards) - dims + 1
    for (r′, r, v) in zip(
        eachslice(advantages, dims = dims),
        eachslice(rewards, dims = dims),
        eachslice(values, dims = dims),
    )
        _generalized_advantage_estimation!(r′, r, v, γ, λ, nothing)
    end
end

function _generalized_advantage_estimation!(
    advantages::AbstractMatrix,
    rewards::AbstractMatrix,
    values::AbstractMatrix,
    γ,
    λ,
    terminal,
    dims::Int,
)
    dims = ndims(rewards) - dims + 1
    for (r′, r, v, t) in zip(
        eachslice(advantages, dims = dims),
        eachslice(rewards, dims = dims),
        eachslice(values, dims = dims),
        eachslice(terminal, dims = dims),
    )
        _generalized_advantage_estimation!(r′, r, v, γ, λ, t)
    end
end

_generalized_advantage_estimation!(
    advantages::AbstractVector,
    rewards::AbstractVector,
    values::AbstractVector,
    γ,
    λ,
    terminal,
    dims::Colon,
) = _generalized_advantage_estimation!(advantages, rewards, values, γ, λ, terminal)

"assuming rewards and advantages are Vector"
function _generalized_advantage_estimation!(advantages, rewards, values, γ, λ, terminal)
    gae = 0
    for i in length(rewards):-1:1
        is_continue = isnothing(terminal) ? true : (!terminal[i])
        delta = rewards[i] + γ * values[i+1] * is_continue - values[i]
        gae = delta + γ * λ * is_continue * gae
        advantages[i] = gae
    end
    advantages
end


text:
module Irradiance
export run_app
using PortAudio, SampledSignals, DSP
import Pkg
const usegpu = "ArrayFire" in keys(Pkg.installed())
using ArrayFire

include("./RotatingBuffer.jl")
include("./LEDTypes.jl")
include("./ParseConfig.jl")
include("./UpdateMethods.jl")
include("./AbstractEffect.jl")
include("./IrradianceCLI.jl")

const old_data = Array{Any, 1}(undef, 0)

function run_app(remote::Bool, args...)
    #ArrayFire.set_backend(AF_BACKEND_OPENCL)
    if remote
        if length(args) >= 1
            args[1]::Int
        end
        if length(args) >= 2
            args[2]::Int
        end
        if length(args) >= 3
            args[3]::IPAddr
        end
        socket = UDPSocket()
        bind(socket, ip"0.0.0.0", 0)
        Sockets.setopt(socket, enable_broadcast=1)
        led_data = remote_config(
            socket,
            length(args)>=1 ? args[1] : 8080,
            length(args)>=2 ? args[2] : 37322,
            length(args)>=3 ? args[3] : ip"10.42.0.255"
        )
    else
        led_data = parse_config(length(args)>=1 && typeof(args[1]) == String ? args[1] : "./lights.json")
        socket = UDPSocket()
    end
    main(led_data, socket)
end

function main(led_data, udpsock)
    signal_channel = Channel{String}(2)
    config_channel = Channel{EffectConfig}(2)
    channels = [signal_channel, config_channel]
    push!(config_channel, EffectConfig(
        HSL(240, 1, 0.5),
        HSL(240, 0, 0.0),
        1,
        1,
        Dict{String, Any}()
    ))
    push!(signal_channel, "1")
    audio = PortAudioStream("default")
    @async handle_cli(channels)
    main_loop(audio, led_data, udpsock, channels)
    close(audio)
    close.(channels)
end

function main_loop(audio, led_data, udpsock, channels)
    signal_channel = channels[1]
    config_channel = channels[2]
    shutdown = false
    clearline = join([" " for i in 1:80])
    for channel in led_data.channels
        channel[1:size(channel,1)] = colorant"black"
    end
    frame_length = (1/10)s
    @sync audio_samp = read(audio, frame_length)
    audio_buf = RotatingBuffer(audio_samp, 2)
    ana = AudioAnalysis(audio_buf[1], 3)
    config = fetch(config_channel)
    current_effect = effect_types["1"](led_data, config, ana)
    current_mode = ""
    while !shutdown
        mode = fetch(signal_channel)
        if fetch(config_channel) != config
            config = fetch(config_channel)
            if issubtype(typeof(current_effect), ConfigurableEffect)
                current_effect.config = config
            end
        end
        if mode != current_mode
            if mode == "shutdown"
                shutdown = true
                current_mode = mode
            elseif mode in keys(effect_types)
                current_effect = effect_types[mode](led_data, config, ana)
                current_mode = mode
            end
        end

        Threads.@threads for i in 1:2
            if Threads.threadid() == 1
                audio_buf[1] = @sync read(audio, frame_length)
            else
                # does the magic of ffts and rotating buffers
                process_audio!(ana, audio_buf[1])

                update!(led_data, ana, current_effect)
            end
        end

        push(led_data, udpsock)
        rotate_clockwise!(audio_buf)
    end
end

function push(led_data::LEDArray, socket::UDPSocket)
    for controller in led_data.controllers
        push(controller, socket)
    end
end

function push(controller::LEDController, socket::UDPSocket)
    controller.raw_data = reshape(controller.addrs', :)
    send(socket, controller.location[1], controller.location[2], controller.raw_data)
end

function convert_to_array(color::ColorTypes.RGB{FixedPointNumbers.Normed{UInt8,8}})
    return [color.r.i, color.g.i, color.b.i]
end

function toHexString(num)
    num <= 255 || error("Input Error")
    magic = "0123456789ABCDEF"
    num1 = floor(Int, num / 16)
    num2 = convert(Int, num % 16)
    return "$(magic[num1+1])$(magic[num2+1])"
end

end


text:
using Test
using DataFrames
using ..ChemistryFeaturization.Atoms: elements
using ..ChemistryFeaturization.FeatureDescriptor: OrbitalFeatureDescriptor
using SparseArrays

@testset "Encode-Decode" begin
    ofd = OrbitalFeatureDescriptor()

    @testset "Encode - with `default_ofd_encode`" begin
        silica = AtomGraph(Float32.([0 1 0; 1 0 1; 0 1 0]), ["O", "Si", "O"])
        encoded_silica = encode(ofd, silica)
        I, J, V = findnz(encoded_silica)

        @test I == [2, 3, 4, 5, 2, 3]
        @test J == [1, 1, 2, 2, 3, 3]
        @test V == [2, 4, 2, 2, 2, 4]

        # non existent element - doesn't really make sense
        nonexistent = AtomGraph(Float32.([1 1; 1 1]), ["Non", "Existent"])
        @test_throws AssertionError ofd(nonexistent)
    end

    @testset "Decode - with `default_ofd_decode`" begin
        alumina = AtomGraph(
            Float32.([0 1 0 0 0; 1 0 1 0 0; 0 1 0 1 0; 0 0 1 0 1; 0 0 0 1 0]),
            ["O", "Al", "O", "Al", "O"],
        )
        @test all(decode(ofd, ofd(alumina)) .== elements(alumina))
    end
end


text:
using WinControllerToFile
using Test

@testset "WinControllerToFile.jl" begin
    subscribe() == nothing
end


text:
using Stochy
using DataStructures

@pp function transition(state)
    state ? flip(0.7) : flip(0.3)
end

@pp function observe(state)
    state ? flip(0.9) : flip(0.1)
end

@pp function hmm(n)
    if n == 0
        local states       = list(true),
              observations = list()
        list(states, observations)
    else
        local prev       = hmm(n-1),
              prevstates = first(prev),
              prevobs    = second(prev),
              newstate   = transition(first(prevstates)),
              newobs     = observe(newstate)
        list(newstate..prevstates, newobs..prevobs)
    end
end

dist = @pp enum() do
    local trueobs = list(false,false,false),
          r = hmm(3)
    factor(trueobs == second(r) ? 0 : -Inf)
    tail(reverse(first(r)))
end

println(dist)


text:
@testset "7.1.2 (d x)^m (a+b arcsinh(c x))^n" begin
    (a, b, c, m, n, x, ) = @variables a b c m n x

    #= ::Package:: =#

    #= ::Title:: =#
    #=Integrands*of*the*form*(d*x)^m*(a+b*asinh(c*x))^n=#

    #= ::Section::Closed:: =#
    #=Integrands*of*the*form*(b*x)^m*asinh(a*x))^n=#

    #= ::Subsection::Closed:: =#
    #=Integrands*of*the*form*x^m*asinh(a*x)^n=#

    #= ::Subsubsection::Closed:: =#
    #=n>0=#

    @test_int [x^4*asinh(a*x), x, 4, -(sqrt(1 + a^2*x^2)/(5*a^5)) + (2*(1 + a^2*x^2)^(3/2))/(15*a^5) - (1 + a^2*x^2)^(5/2)/(25*a^5) + (1/5)*x^5*asinh(a*x)]
    @test_int [x^3*asinh(a*x), x, 4, (3*x*sqrt(1 + a^2*x^2))/(32*a^3) - (x^3*sqrt(1 + a^2*x^2))/(16*a) - (3*asinh(a*x))/(32*a^4) + (1/4)*x^4*asinh(a*x)]
    @test_int [x^2*asinh(a*x), x, 4, sqrt(1 + a^2*x^2)/(3*a^3) - (1 + a^2*x^2)^(3/2)/(9*a^3) + (1/3)*x^3*asinh(a*x)]
    @test_int [x^1*asinh(a*x), x, 3, -((x*sqrt(1 + a^2*x^2))/(4*a)) + asinh(a*x)/(4*a^2) + (1/2)*x^2*asinh(a*x)]
    @test_int [x^0*asinh(a*x), x, 2, -(sqrt(1 + a^2*x^2)/a) + x*asinh(a*x)]
    @test_int [asinh(a*x)/x^1, x, 5, (-(1/2))*asinh(a*x)^2 + asinh(a*x)*log(1 - ℯ^(2*asinh(a*x))) + (1/2)*PolyLog(2, ℯ^(2*asinh(a*x)))]
    @test_int [asinh(a*x)/x^2, x, 4, -(asinh(a*x)/x) - a*atanh(sqrt(1 + a^2*x^2))]
    @test_int [asinh(a*x)/x^3, x, 2, -((a*sqrt(1 + a^2*x^2))/(2*x)) - asinh(a*x)/(2*x^2)]
    @test_int [asinh(a*x)/x^4, x, 5, -((a*sqrt(1 + a^2*x^2))/(6*x^2)) - asinh(a*x)/(3*x^3) + (1/6)*a^3*atanh(sqrt(1 + a^2*x^2))]
    @test_int [asinh(a*x)/x^5, x, 3, -((a*sqrt(1 + a^2*x^2))/(12*x^3)) + (a^3*sqrt(1 + a^2*x^2))/(6*x) - asinh(a*x)/(4*x^4)]
    @test_int [asinh(a*x)/x^6, x, 6, -((a*sqrt(1 + a^2*x^2))/(20*x^4)) + (3*a^3*sqrt(1 + a^2*x^2))/(40*x^2) - asinh(a*x)/(5*x^5) - (3/40)*a^5*atanh(sqrt(1 + a^2*x^2))]

    @test_int [x^4*asinh(a*x)^2, x, 7, (16*x)/(75*a^4) - (8*x^3)/(225*a^2) + (2*x^5)/125 - (16*sqrt(1 + a^2*x^2)*asinh(a*x))/(75*a^5) + (8*x^2*sqrt(1 + a^2*x^2)*asinh(a*x))/(75*a^3) - (2*x^4*sqrt(1 + a^2*x^2)*asinh(a*x))/(25*a) + (1/5)*x^5*asinh(a*x)^2]
    @test_int [x^3*asinh(a*x)^2, x, 6, (-3*x^2)/(32*a^2) + x^4/32 + (3*x*sqrt(1 + a^2*x^2)*asinh(a*x))/(16*a^3) - (x^3*sqrt(1 + a^2*x^2)*asinh(a*x))/(8*a) - (3*asinh(a*x)^2)/(32*a^4) + (x^4*asinh(a*x)^2)/4]
    @test_int [x^2*asinh(a*x)^2, x, 5, -((4*x)/(9*a^2)) + (2*x^3)/27 + (4*sqrt(1 + a^2*x^2)*asinh(a*x))/(9*a^3) - (2*x^2*sqrt(1 + a^2*x^2)*asinh(a*x))/(9*a) + (1/3)*x^3*asinh(a*x)^2]
    @test_int [x*asinh(a*x)^2, x, 4, x^2/4 - (x*sqrt(1 + a^2*x^2)*asinh(a*x))/(2*a) + asinh(a*x)^2/(4*a^2) + (x^2*asinh(a*x)^2)/2]
    @test_int [asinh(a*x)^2, x, 3, 2*x - (2*sqrt(1 + a^2*x^2)*asinh(a*x))/a + x*asinh(a*x)^2]
    @test_int [asinh(a*x)^2/x, x, 6, -asinh(a*x)^3/3 + asinh(a*x)^2*log(1 - ℯ^(2*asinh(a*x))) + asinh(a*x)*PolyLog(2, ℯ^(2*asinh(a*x))) - PolyLog(3, ℯ^(2*asinh(a*x)))/2]
    @test_int [asinh(a*x)^2/x^2, x, 7, -(asinh(a*x)^2/x) - 4*a*asinh(a*x)*atanh(ℯ^asinh(a*x)) - 2*a*PolyLog(2, -ℯ^asinh(a*x)) + 2*a*PolyLog(2, ℯ^asinh(a*x))]
    @test_int [asinh(a*x)^2/x^3, x, 3, -((a*sqrt(1 + a^2*x^2)*asinh(a*x))/x) - asinh(a*x)^2/(2*x^2) + a^2*log(x)]
    @test_int [asinh(a*x)^2/x^4, x, 9, -(a^2/(3*x)) - (a*sqrt(1 + a^2*x^2)*asinh(a*x))/(3*x^2) - asinh(a*x)^2/(3*x^3) + (2/3)*a^3*asinh(a*x)*atanh(ℯ^asinh(a*x)) + (1/3)*a^3*PolyLog(2, -ℯ^asinh(a*x)) - (1/3)*a^3*PolyLog(2, ℯ^asinh(a*x))]
    @test_int [asinh(a*x)^2/x^5, x, 5, -a^2/(12*x^2) - (a*sqrt(1 + a^2*x^2)*asinh(a*x))/(6*x^3) + (a^3*sqrt(1 + a^2*x^2)*asinh(a*x))/(3*x) - asinh(a*x)^2/(4*x^4) - (a^4*log(x))/3]

    @test_int [x^4*asinh(a*x)^3, x, 14, -((298*sqrt(1 + a^2*x^2))/(375*a^5)) + (76*(1 + a^2*x^2)^(3/2))/(1125*a^5) - (6*(1 + a^2*x^2)^(5/2))/(625*a^5) + (16*x*asinh(a*x))/(25*a^4) - (8*x^3*asinh(a*x))/(75*a^2) + (6/125)*x^5*asinh(a*x) - (8*sqrt(1 + a^2*x^2)*asinh(a*x)^2)/(25*a^5) + (4*x^2*sqrt(1 + a^2*x^2)*asinh(a*x)^2)/(25*a^3) - (3*x^4*sqrt(1 + a^2*x^2)*asinh(a*x)^2)/(25*a) + (1/5)*x^5*asinh(a*x)^3]
    @test_int [x^3*asinh(a*x)^3, x, 11, (45*x*sqrt(1 + a^2*x^2))/(256*a^3) - (3*x^3*sqrt(1 + a^2*x^2))/(128*a) - (45*asinh(a*x))/(256*a^4) - (9*x^2*asinh(a*x))/(32*a^2) + (3*x^4*asinh(a*x))/32 + (9*x*sqrt(1 + a^2*x^2)*asinh(a*x)^2)/(32*a^3) - (3*x^3*sqrt(1 + a^2*x^2)*asinh(a*x)^2)/(16*a) - (3*asinh(a*x)^3)/(32*a^4) + (x^4*asinh(a*x)^3)/4]
    @test_int [x^2*asinh(a*x)^3, x, 9, (14*sqrt(1 + a^2*x^2))/(9*a^3) - (2*(1 + a^2*x^2)^(3/2))/(27*a^3) - (4*x*asinh(a*x))/(3*a^2) + (2/9)*x^3*asinh(a*x) + (2*sqrt(1 + a^2*x^2)*asinh(a*x)^2)/(3*a^3) - (x^2*sqrt(1 + a^2*x^2)*asinh(a*x)^2)/(3*a) + (1/3)*x^3*asinh(a*x)^3]
    @test_int [x*asinh(a*x)^3, x, 6, (-3*x*sqrt(1 + a^2*x^2))/(8*a) + (3*asinh(a*x))/(8*a^2) + (3*x^2*asinh(a*x))/4 - (3*x*sqrt(1 + a^2*x^2)*asinh(a*x)^2)/(4*a) + asinh(a*x)^3/(4*a^2) + (x^2*asinh(a*x)^3)/2]
    @test_int [asinh(a*x)^3, x, 4, (-6*sqrt(1 + a^2*x^2))/a + 6*x*asinh(a*x) - (3*sqrt(1 + a^2*x^2)*asinh(a*x)^2)/a + x*asinh(a*x)^3]
    @test_int [asinh(a*x)^3/x, x, 7, -asinh(a*x)^4/4 + asinh(a*x)^3*log(1 - ℯ^(2*asinh(a*x))) + (3*asinh(a*x)^2*PolyLog(2, ℯ^(2*asinh(a*x))))/2 - (3*asinh(a*x)*PolyLog(3, ℯ^(2*asinh(a*x))))/2 + (3*PolyLog(4, ℯ^(2*asinh(a*x))))/4]
    @test_int [asinh(a*x)^3/x^2, x, 9, -(asinh(a*x)^3/x) - 6*a*asinh(a*x)^2*atanh(ℯ^asinh(a*x)) - 6*a*asinh(a*x)*PolyLog(2, -ℯ^asinh(a*x)) + 6*a*asinh(a*x)*PolyLog(2, ℯ^asinh(a*x)) + 6*a*PolyLog(3, -ℯ^asinh(a*x)) - 6*a*PolyLog(3, ℯ^asinh(a*x))]
    @test_int [asinh(a*x)^3/x^3, x, 7, (-3*a^2*asinh(a*x)^2)/2 - (3*a*sqrt(1 + a^2*x^2)*asinh(a*x)^2)/(2*x) - asinh(a*x)^3/(2*x^2) + 3*a^2*asinh(a*x)*log(1 - ℯ^(2*asinh(a*x))) + (3*a^2*PolyLog(2, ℯ^(2*asinh(a*x))))/2]
    @test_int [asinh(a*x)^3/x^4, x, 14, -((a^2*asinh(a*x))/x) - (a*sqrt(1 + a^2*x^2)*asinh(a*x)^2)/(2*x^2) - asinh(a*x)^3/(3*x^3) + a^3*asinh(a*x)^2*atanh(ℯ^asinh(a*x)) - a^3*atanh(sqrt(1 + a^2*x^2)) + a^3*asinh(a*x)*PolyLog(2, -ℯ^asinh(a*x)) - a^3*asinh(a*x)*PolyLog(2, ℯ^asinh(a*x)) - a^3*PolyLog(3, -ℯ^asinh(a*x)) + a^3*PolyLog(3, ℯ^asinh(a*x))]
    @test_int [asinh(a*x)^3/x^5, x, 10, -(a^3*sqrt(1 + a^2*x^2))/(4*x) - (a^2*asinh(a*x))/(4*x^2) + (a^4*asinh(a*x)^2)/2 - (a*sqrt(1 + a^2*x^2)*asinh(a*x)^2)/(4*x^3) + (a^3*sqrt(1 + a^2*x^2)*asinh(a*x)^2)/(2*x) - asinh(a*x)^3/(4*x^4) - a^4*asinh(a*x)*log(1 - ℯ^(2*asinh(a*x))) - (a^4*PolyLog(2, ℯ^(2*asinh(a*x))))/2]

    @test_int [x^5*asinh(a*x)^4, x, 23, (245*x^2)/(1152*a^4) - (65*x^4)/(3456*a^2) + x^6/324 - (245*x*sqrt(1 + a^2*x^2)*asinh(a*x))/(576*a^5) + (65*x^3*sqrt(1 + a^2*x^2)*asinh(a*x))/(864*a^3) - (x^5*sqrt(1 + a^2*x^2)*asinh(a*x))/(54*a) + (245*asinh(a*x)^2)/(1152*a^6) + (5*x^2*asinh(a*x)^2)/(16*a^4) - (5*x^4*asinh(a*x)^2)/(48*a^2) + (x^6*asinh(a*x)^2)/18 - (5*x*sqrt(1 + a^2*x^2)*asinh(a*x)^3)/(24*a^5) + (5*x^3*sqrt(1 + a^2*x^2)*asinh(a*x)^3)/(36*a^3) - (x^5*sqrt(1 + a^2*x^2)*asinh(a*x)^3)/(9*a) + (5*asinh(a*x)^4)/(96*a^6) + (x^6*asinh(a*x)^4)/6]
    @test_int [x^4*asinh(a*x)^4, x, 19, (16576*x)/(5625*a^4) - (1088*x^3)/(16875*a^2) + (24*x^5)/3125 - (16576*sqrt(1 + a^2*x^2)*asinh(a*x))/(5625*a^5) + (1088*x^2*sqrt(1 + a^2*x^2)*asinh(a*x))/(5625*a^3) - (24*x^4*sqrt(1 + a^2*x^2)*asinh(a*x))/(625*a) + (32*x*asinh(a*x)^2)/(25*a^4) - (16*x^3*asinh(a*x)^2)/(75*a^2) + (12/125)*x^5*asinh(a*x)^2 - (32*sqrt(1 + a^2*x^2)*asinh(a*x)^3)/(75*a^5) + (16*x^2*sqrt(1 + a^2*x^2)*asinh(a*x)^3)/(75*a^3) - (4*x^4*sqrt(1 + a^2*x^2)*asinh(a*x)^3)/(25*a) + (1/5)*x^5*asinh(a*x)^4]
    @test_int [x^3*asinh(a*x)^4, x, 14, (-45*x^2)/(128*a^2) + (3*x^4)/128 + (45*x*sqrt(1 + a^2*x^2)*asinh(a*x))/(64*a^3) - (3*x^3*sqrt(1 + a^2*x^2)*asinh(a*x))/(32*a) - (45*asinh(a*x)^2)/(128*a^4) - (9*x^2*asinh(a*x)^2)/(16*a^2) + (3*x^4*asinh(a*x)^2)/16 + (3*x*sqrt(1 + a^2*x^2)*asinh(a*x)^3)/(8*a^3) - (x^3*sqrt(1 + a^2*x^2)*asinh(a*x)^3)/(4*a) - (3*asinh(a*x)^4)/(32*a^4) + (x^4*asinh(a*x)^4)/4]
    @test_int [x^2*asinh(a*x)^4, x, 11, -((160*x)/(27*a^2)) + (8*x^3)/81 + (160*sqrt(1 + a^2*x^2)*asinh(a*x))/(27*a^3) - (8*x^2*sqrt(1 + a^2*x^2)*asinh(a*x))/(27*a) - (8*x*asinh(a*x)^2)/(3*a^2) + (4/9)*x^3*asinh(a*x)^2 + (8*sqrt(1 + a^2*x^2)*asinh(a*x)^3)/(9*a^3) - (4*x^2*sqrt(1 + a^2*x^2)*asinh(a*x)^3)/(9*a) + (1/3)*x^3*asinh(a*x)^4]
    @test_int [x*asinh(a*x)^4, x, 7, (3*x^2)/4 - (3*x*sqrt(1 + a^2*x^2)*asinh(a*x))/(2*a) + (3*asinh(a*x)^2)/(4*a^2) + (3*x^2*asinh(a*x)^2)/2 - (x*sqrt(1 + a^2*x^2)*asinh(a*x)^3)/a + asinh(a*x)^4/(4*a^2) + (x^2*asinh(a*x)^4)/2]
    @test_int [asinh(a*x)^4, x, 5, 24*x - (24*sqrt(1 + a^2*x^2)*asinh(a*x))/a + 12*x*asinh(a*x)^2 - (4*sqrt(1 + a^2*x^2)*asinh(a*x)^3)/a + x*asinh(a*x)^4]
    @test_int [asinh(a*x)^4/x, x, 8, -asinh(a*x)^5/5 + asinh(a*x)^4*log(1 - ℯ^(2*asinh(a*x))) + 2*asinh(a*x)^3*PolyLog(2, ℯ^(2*asinh(a*x))) - 3*asinh(a*x)^2*PolyLog(3, ℯ^(2*asinh(a*x))) + 3*asinh(a*x)*PolyLog(4, ℯ^(2*asinh(a*x))) - (3*PolyLog(5, ℯ^(2*asinh(a*x))))/2]
    @test_int [asinh(a*x)^4/x^2, x, 11, -(asinh(a*x)^4/x) - 8*a*asinh(a*x)^3*atanh(ℯ^asinh(a*x)) - 12*a*asinh(a*x)^2*PolyLog(2, -ℯ^asinh(a*x)) + 12*a*asinh(a*x)^2*PolyLog(2, ℯ^asinh(a*x)) + 24*a*asinh(a*x)*PolyLog(3, -ℯ^asinh(a*x)) - 24*a*asinh(a*x)*PolyLog(3, ℯ^asinh(a*x)) - 24*a*PolyLog(4, -ℯ^asinh(a*x)) + 24*a*PolyLog(4, ℯ^asinh(a*x))]
    @test_int [asinh(a*x)^4/x^3, x, 8, -2*a^2*asinh(a*x)^3 - (2*a*sqrt(1 + a^2*x^2)*asinh(a*x)^3)/x - asinh(a*x)^4/(2*x^2) + 6*a^2*asinh(a*x)^2*log(1 - ℯ^(2*asinh(a*x))) + 6*a^2*asinh(a*x)*PolyLog(2, ℯ^(2*asinh(a*x))) - 3*a^2*PolyLog(3, ℯ^(2*asinh(a*x)))]
    @test_int [asinh(a*x)^4/x^4, x, 19, -((2*a^2*asinh(a*x)^2)/x) - (2*a*sqrt(1 + a^2*x^2)*asinh(a*x)^3)/(3*x^2) - asinh(a*x)^4/(3*x^3) - 8*a^3*asinh(a*x)*atanh(ℯ^asinh(a*x)) + (4/3)*a^3*asinh(a*x)^3*atanh(ℯ^asinh(a*x)) - 4*a^3*PolyLog(2, -ℯ^asinh(a*x)) + 2*a^3*asinh(a*x)^2*PolyLog(2, -ℯ^asinh(a*x)) + 4*a^3*PolyLog(2, ℯ^asinh(a*x)) - 2*a^3*asinh(a*x)^2*PolyLog(2, ℯ^asinh(a*x)) - 4*a^3*asinh(a*x)*PolyLog(3, -ℯ^asinh(a*x)) + 4*a^3*asinh(a*x)*PolyLog(3, ℯ^asinh(a*x)) + 4*a^3*PolyLog(4, -ℯ^asinh(a*x)) - 4*a^3*PolyLog(4, ℯ^asinh(a*x))]

    #= ::Subsubsection::Closed:: =#
    #=n<0=#

    @test_int [x^6/asinh(a*x), x, 7, -((5*CoshIntegral(asinh(a*x)))/(64*a^7)) + (9*CoshIntegral(3*asinh(a*x)))/(64*a^7) - (5*CoshIntegral(5*asinh(a*x)))/(64*a^7) + CoshIntegral(7*asinh(a*x))/(64*a^7)]
    @test_int [x^5/asinh(a*x), x, 6, (5*SinhIntegral(2*asinh(a*x)))/(32*a^6) - SinhIntegral(4*asinh(a*x))/(8*a^6) + SinhIntegral(6*asinh(a*x))/(32*a^6)]
    @test_int [x^4/asinh(a*x), x, 6, CoshIntegral(asinh(a*x))/(8*a^5) - (3*CoshIntegral(3*asinh(a*x)))/(16*a^5) + CoshIntegral(5*asinh(a*x))/(16*a^5)]
    @test_int [x^3/asinh(a*x), x, 5, -SinhIntegral(2*asinh(a*x))/(4*a^4) + SinhIntegral(4*asinh(a*x))/(8*a^4)]
    @test_int [x^2/asinh(a*x), x, 5, -CoshIntegral(asinh(a*x))/(4*a^3) + CoshIntegral(3*asinh(a*x))/(4*a^3)]
    @test_int [x/asinh(a*x), x, 4, SinhIntegral(2*asinh(a*x))/(2*a^2)]
    @test_int [asinh(a*x)^(-1), x, 2, CoshIntegral(asinh(a*x))/a]
    @test_int [1/(x*asinh(a*x)), x, 0, Unintegrable(1/(x*asinh(a*x)), x)]
    @test_int [1/(x^2*asinh(a*x)), x, 0, Unintegrable(1/(x^2*asinh(a*x)), x)]

    @test_int [x^6/asinh(a*x)^2, x, 6, -((x^6*sqrt(1 + a^2*x^2))/(a*asinh(a*x))) - (5*SinhIntegral(asinh(a*x)))/(64*a^7) + (27*SinhIntegral(3*asinh(a*x)))/(64*a^7) - (25*SinhIntegral(5*asinh(a*x)))/(64*a^7) + (7*SinhIntegral(7*asinh(a*x)))/(64*a^7)]
    @test_int [x^5/asinh(a*x)^2, x, 5, -((x^5*sqrt(1 + a^2*x^2))/(a*asinh(a*x))) + (5*CoshIntegral(2*asinh(a*x)))/(16*a^6) - CoshIntegral(4*asinh(a*x))/(2*a^6) + (3*CoshIntegral(6*asinh(a*x)))/(16*a^6)]
    @test_int [x^4/asinh(a*x)^2, x, 5, -((x^4*sqrt(1 + a^2*x^2))/(a*asinh(a*x))) + SinhIntegral(asinh(a*x))/(8*a^5) - (9*SinhIntegral(3*asinh(a*x)))/(16*a^5) + (5*SinhIntegral(5*asinh(a*x)))/(16*a^5)]
    @test_int [x^3/asinh(a*x)^2, x, 4, -((x^3*sqrt(1 + a^2*x^2))/(a*asinh(a*x))) - CoshIntegral(2*asinh(a*x))/(2*a^4) + CoshIntegral(4*asinh(a*x))/(2*a^4)]
    @test_int [x^2/asinh(a*x)^2, x, 4, -((x^2*sqrt(1 + a^2*x^2))/(a*asinh(a*x))) - SinhIntegral(asinh(a*x))/(4*a^3) + (3*SinhIntegral(3*asinh(a*x)))/(4*a^3)]
    @test_int [x/asinh(a*x)^2, x, 2, -((x*sqrt(1 + a^2*x^2))/(a*asinh(a*x))) + CoshIntegral(2*asinh(a*x))/a^2]
    @test_int [asinh(a*x)^(-2), x, 3, -(sqrt(1 + a^2*x^2)/(a*asinh(a*x))) + SinhIntegral(asinh(a*x))/a]
    @test_int [1/(x*asinh(a*x)^2), x, 0, Unintegrable(1/(x*asinh(a*x)^2), x)]
    @test_int [1/(x^2*asinh(a*x)^2), x, 0, Unintegrable(1/(x^2*asinh(a*x)^2), x)]

    @test_int [x^4/asinh(a*x)^3, x, 14, -(x^4*sqrt(1 + a^2*x^2))/(2*a*asinh(a*x)^2) - (2*x^3)/(a^2*asinh(a*x)) - (5*x^5)/(2*asinh(a*x)) + CoshIntegral(asinh(a*x))/(16*a^5) - (27*CoshIntegral(3*asinh(a*x)))/(32*a^5) + (25*CoshIntegral(5*asinh(a*x)))/(32*a^5)]
    @test_int [x^3/asinh(a*x)^3, x, 12, -(x^3*sqrt(1 + a^2*x^2))/(2*a*asinh(a*x)^2) - (3*x^2)/(2*a^2*asinh(a*x)) - (2*x^4)/asinh(a*x) - SinhIntegral(2*asinh(a*x))/(2*a^4) + SinhIntegral(4*asinh(a*x))/a^4]
    @test_int [x^2/asinh(a*x)^3, x, 10, -(x^2*sqrt(1 + a^2*x^2))/(2*a*asinh(a*x)^2) - x/(a^2*asinh(a*x)) - (3*x^3)/(2*asinh(a*x)) - CoshIntegral(asinh(a*x))/(8*a^3) + (9*CoshIntegral(3*asinh(a*x)))/(8*a^3)]
    @test_int [x/asinh(a*x)^3, x, 7, -(x*sqrt(1 + a^2*x^2))/(2*a*asinh(a*x)^2) - 1/(2*a^2*asinh(a*x)) - x^2/asinh(a*x) + SinhIntegral(2*asinh(a*x))/a^2]
    @test_int [asinh(a*x)^(-3), x, 4, -sqrt(1 + a^2*x^2)/(2*a*asinh(a*x)^2) - x/(2*asinh(a*x)) + CoshIntegral(asinh(a*x))/(2*a)]
    @test_int [1/(x*asinh(a*x)^3), x, 0, Unintegrable(1/(x*asinh(a*x)^3), x)]
    @test_int [1/(x^2*asinh(a*x)^3), x, 0, Unintegrable(1/(x^2*asinh(a*x)^3), x)]

    @test_int [x^4/asinh(a*x)^4, x, 12, -(x^4*sqrt(1 + a^2*x^2))/(3*a*asinh(a*x)^3) - (2*x^3)/(3*a^2*asinh(a*x)^2) - (5*x^5)/(6*asinh(a*x)^2) - (2*x^2*sqrt(1 + a^2*x^2))/(a^3*asinh(a*x)) - (25*x^4*sqrt(1 + a^2*x^2))/(6*a*asinh(a*x)) + SinhIntegral(asinh(a*x))/(48*a^5) - (27*SinhIntegral(3*asinh(a*x)))/(32*a^5) + (125*SinhIntegral(5*asinh(a*x)))/(96*a^5)]
    @test_int [x^3/asinh(a*x)^4, x, 9, -(x^3*sqrt(1 + a^2*x^2))/(3*a*asinh(a*x)^3) - x^2/(2*a^2*asinh(a*x)^2) - (2*x^4)/(3*asinh(a*x)^2) - (x*sqrt(1 + a^2*x^2))/(a^3*asinh(a*x)) - (8*x^3*sqrt(1 + a^2*x^2))/(3*a*asinh(a*x)) - CoshIntegral(2*asinh(a*x))/(3*a^4) + (4*CoshIntegral(4*asinh(a*x)))/(3*a^4)]
    @test_int [x^2/asinh(a*x)^4, x, 10, -(x^2*sqrt(1 + a^2*x^2))/(3*a*asinh(a*x)^3) - x/(3*a^2*asinh(a*x)^2) - x^3/(2*asinh(a*x)^2) - sqrt(1 + a^2*x^2)/(3*a^3*asinh(a*x)) - (3*x^2*sqrt(1 + a^2*x^2))/(2*a*asinh(a*x)) - SinhIntegral(asinh(a*x))/(24*a^3) + (9*SinhIntegral(3*asinh(a*x)))/(8*a^3)]
    @test_int [x/asinh(a*x)^4, x, 5, -(x*sqrt(1 + a^2*x^2))/(3*a*asinh(a*x)^3) - 1/(6*a^2*asinh(a*x)^2) - x^2/(3*asinh(a*x)^2) - (2*x*sqrt(1 + a^2*x^2))/(3*a*asinh(a*x)) + (2*CoshIntegral(2*asinh(a*x)))/(3*a^2)]
    @test_int [asinh(a*x)^(-4), x, 5, -sqrt(1 + a^2*x^2)/(3*a*asinh(a*x)^3) - x/(6*asinh(a*x)^2) - sqrt(1 + a^2*x^2)/(6*a*asinh(a*x)) + SinhIntegral(asinh(a*x))/(6*a)]
    @test_int [1/(x*asinh(a*x)^4), x, 0, Unintegrable(1/(x*asinh(a*x)^4), x)]
    @test_int [1/(x^2*asinh(a*x)^4), x, 0, Unintegrable(1/(x^2*asinh(a*x)^4), x)]

    #= ::Subsection::Closed:: =#
    #=Integrands*of*the*form*x^m*asinh(a*x)^(n/2)=#

    #= ::Subsubsection::Closed:: =#
    #=n>0=#

    @test_int [x^4*sqrt(asinh(a*x)), x, 19, (x^5*sqrt(asinh(a*x)))/5 + (sqrt(pi)*erf(sqrt(asinh(a*x))))/(32*a^5) - (sqrt(pi/3)*erf(sqrt(3)*sqrt(asinh(a*x))))/(64*a^5) + (sqrt(pi/5)*erf(sqrt(5)*sqrt(asinh(a*x))))/(320*a^5) - (sqrt(pi)*erfi(sqrt(asinh(a*x))))/(32*a^5) + (sqrt(pi/3)*erfi(sqrt(3)*sqrt(asinh(a*x))))/(64*a^5) - (sqrt(pi/5)*erfi(sqrt(5)*sqrt(asinh(a*x))))/(320*a^5)]
    @test_int [x^3*sqrt(asinh(a*x)), x, 14, (-3*sqrt(asinh(a*x)))/(32*a^4) + (x^4*sqrt(asinh(a*x)))/4 - (sqrt(pi)*erf(2*sqrt(asinh(a*x))))/(256*a^4) + (sqrt(pi/2)*erf(sqrt(2)*sqrt(asinh(a*x))))/(32*a^4) - (sqrt(pi)*erfi(2*sqrt(asinh(a*x))))/(256*a^4) + (sqrt(pi/2)*erfi(sqrt(2)*sqrt(asinh(a*x))))/(32*a^4)]
    @test_int [x^2*sqrt(asinh(a*x)), x, 14, (x^3*sqrt(asinh(a*x)))/3 - (sqrt(pi)*erf(sqrt(asinh(a*x))))/(16*a^3) + (sqrt(pi/3)*erf(sqrt(3)*sqrt(asinh(a*x))))/(48*a^3) + (sqrt(pi)*erfi(sqrt(asinh(a*x))))/(16*a^3) - (sqrt(pi/3)*erfi(sqrt(3)*sqrt(asinh(a*x))))/(48*a^3)]
    @test_int [x*sqrt(asinh(a*x)), x, 9, sqrt(asinh(a*x))/(4*a^2) + (x^2*sqrt(asinh(a*x)))/2 - (sqrt(pi/2)*erf(sqrt(2)*sqrt(asinh(a*x))))/(16*a^2) - (sqrt(pi/2)*erfi(sqrt(2)*sqrt(asinh(a*x))))/(16*a^2)]
    @test_int [sqrt(asinh(a*x)), x, 7, x*sqrt(asinh(a*x)) + (sqrt(pi)*erf(sqrt(asinh(a*x))))/(4*a) - (sqrt(pi)*erfi(sqrt(asinh(a*x))))/(4*a)]
    @test_int [sqrt(asinh(a*x))/x, x, 0, Unintegrable(sqrt(asinh(a*x))/x, x)]

    @test_int [x^4*asinh(a*x)^(3/2), x, 41, -((4*sqrt(1 + a^2*x^2)*sqrt(asinh(a*x)))/(25*a^5)) + (2*x^2*sqrt(1 + a^2*x^2)*sqrt(asinh(a*x)))/(25*a^3) - (3*x^4*sqrt(1 + a^2*x^2)*sqrt(asinh(a*x)))/(50*a) + (1/5)*x^5*asinh(a*x)^(3/2) + (3*sqrt(pi)*erf(sqrt(asinh(a*x))))/(64*a^5) - (sqrt(pi/3)*erf(sqrt(3)*sqrt(asinh(a*x))))/(200*a^5) - (3*sqrt(3*pi)*erf(sqrt(3)*sqrt(asinh(a*x))))/(3200*a^5) + (3*sqrt(pi/5)*erf(sqrt(5)*sqrt(asinh(a*x))))/(3200*a^5) + (3*sqrt(pi)*erfi(sqrt(asinh(a*x))))/(64*a^5) - (sqrt(pi/3)*erfi(sqrt(3)*sqrt(asinh(a*x))))/(200*a^5) - (3*sqrt(3*pi)*erfi(sqrt(3)*sqrt(asinh(a*x))))/(3200*a^5) + (3*sqrt(pi/5)*erfi(sqrt(5)*sqrt(asinh(a*x))))/(3200*a^5)]
    @test_int [x^3*asinh(a*x)^(3/2), x, 25, (9*x*sqrt(1 + a^2*x^2)*sqrt(asinh(a*x)))/(64*a^3) - (3*x^3*sqrt(1 + a^2*x^2)*sqrt(asinh(a*x)))/(32*a) - (3*asinh(a*x)^(3/2))/(32*a^4) + (1/4)*x^4*asinh(a*x)^(3/2) - (3*sqrt(pi)*erf(2*sqrt(asinh(a*x))))/(2048*a^4) + (3*sqrt(pi/2)*erf(sqrt(2)*sqrt(asinh(a*x))))/(128*a^4) + (3*sqrt(pi)*erfi(2*sqrt(asinh(a*x))))/(2048*a^4) - (3*sqrt(pi/2)*erfi(sqrt(2)*sqrt(asinh(a*x))))/(128*a^4)]
    @test_int [x^2*asinh(a*x)^(3/2), x, 22, (sqrt(1 + a^2*x^2)*sqrt(asinh(a*x)))/(3*a^3) - (x^2*sqrt(1 + a^2*x^2)*sqrt(asinh(a*x)))/(6*a) + (1/3)*x^3*asinh(a*x)^(3/2) - (3*sqrt(pi)*erf(sqrt(asinh(a*x))))/(32*a^3) + (sqrt(pi/3)*erf(sqrt(3)*sqrt(asinh(a*x))))/(96*a^3) - (3*sqrt(pi)*erfi(sqrt(asinh(a*x))))/(32*a^3) + (sqrt(pi/3)*erfi(sqrt(3)*sqrt(asinh(a*x))))/(96*a^3)]
    @test_int [x^1*asinh(a*x)^(3/2), x, 11, -((3*x*sqrt(1 + a^2*x^2)*sqrt(asinh(a*x)))/(8*a)) + asinh(a*x)^(3/2)/(4*a^2) + (1/2)*x^2*asinh(a*x)^(3/2) - (3*sqrt(pi/2)*erf(sqrt(2)*sqrt(asinh(a*x))))/(64*a^2) + (3*sqrt(pi/2)*erfi(sqrt(2)*sqrt(asinh(a*x))))/(64*a^2)]
    @test_int [asinh(a*x)^(3/2), x, 8, (-3*sqrt(1 + a^2*x^2)*sqrt(asinh(a*x)))/(2*a) + x*asinh(a*x)^(3/2) + (3*sqrt(pi)*erf(sqrt(asinh(a*x))))/(8*a) + (3*sqrt(pi)*erfi(sqrt(asinh(a*x))))/(8*a)]
    @test_int [asinh(a*x)^(3/2)/x, x, 0, Unintegrable(asinh(a*x)^(3/2)/x, x)]

    @test_int [x^4*asinh(a*x)^(5/2), x, 44, (2*x*sqrt(asinh(a*x)))/(5*a^4) - (x^3*sqrt(asinh(a*x)))/(15*a^2) + (3/100)*x^5*sqrt(asinh(a*x)) - (4*sqrt(1 + a^2*x^2)*asinh(a*x)^(3/2))/(15*a^5) + (2*x^2*sqrt(1 + a^2*x^2)*asinh(a*x)^(3/2))/(15*a^3) - (x^4*sqrt(1 + a^2*x^2)*asinh(a*x)^(3/2))/(10*a) + (1/5)*x^5*asinh(a*x)^(5/2) + (15*sqrt(pi)*erf(sqrt(asinh(a*x))))/(128*a^5) - (sqrt(pi/3)*erf(sqrt(3)*sqrt(asinh(a*x))))/(240*a^5) - (sqrt(3*pi)*erf(sqrt(3)*sqrt(asinh(a*x))))/(1280*a^5) + (3*sqrt(pi/5)*erf(sqrt(5)*sqrt(asinh(a*x))))/(6400*a^5) - (15*sqrt(pi)*erfi(sqrt(asinh(a*x))))/(128*a^5) + (sqrt(pi/3)*erfi(sqrt(3)*sqrt(asinh(a*x))))/(240*a^5) + (sqrt(3*pi)*erfi(sqrt(3)*sqrt(asinh(a*x))))/(1280*a^5) - (3*sqrt(pi/5)*erfi(sqrt(5)*sqrt(asinh(a*x))))/(6400*a^5)]
    @test_int [x^3*asinh(a*x)^(5/2), x, 27, (-225*sqrt(asinh(a*x)))/(2048*a^4) - (45*x^2*sqrt(asinh(a*x)))/(256*a^2) + (15*x^4*sqrt(asinh(a*x)))/256 + (15*x*sqrt(1 + a^2*x^2)*asinh(a*x)^(3/2))/(64*a^3) - (5*x^3*sqrt(1 + a^2*x^2)*asinh(a*x)^(3/2))/(32*a) - (3*asinh(a*x)^(5/2))/(32*a^4) + (x^4*asinh(a*x)^(5/2))/4 - (15*sqrt(pi)*erf(2*sqrt(asinh(a*x))))/(16384*a^4) + (15*sqrt(pi/2)*erf(sqrt(2)*sqrt(asinh(a*x))))/(512*a^4) - (15*sqrt(pi)*erfi(2*sqrt(asinh(a*x))))/(16384*a^4) + (15*sqrt(pi/2)*erfi(sqrt(2)*sqrt(asinh(a*x))))/(512*a^4)]
    @test_int [x^2*asinh(a*x)^(5/2), x, 24, (-5*x*sqrt(asinh(a*x)))/(6*a^2) + (5*x^3*sqrt(asinh(a*x)))/36 + (5*sqrt(1 + a^2*x^2)*asinh(a*x)^(3/2))/(9*a^3) - (5*x^2*sqrt(1 + a^2*x^2)*asinh(a*x)^(3/2))/(18*a) + (x^3*asinh(a*x)^(5/2))/3 - (15*sqrt(pi)*erf(sqrt(asinh(a*x))))/(64*a^3) + (5*sqrt(pi/3)*erf(sqrt(3)*sqrt(asinh(a*x))))/(576*a^3) + (15*sqrt(pi)*erfi(sqrt(asinh(a*x))))/(64*a^3) - (5*sqrt(pi/3)*erfi(sqrt(3)*sqrt(asinh(a*x))))/(576*a^3)]
    @test_int [x*asinh(a*x)^(5/2), x, 12, (15*sqrt(asinh(a*x)))/(64*a^2) + (15*x^2*sqrt(asinh(a*x)))/32 - (5*x*sqrt(1 + a^2*x^2)*asinh(a*x)^(3/2))/(8*a) + asinh(a*x)^(5/2)/(4*a^2) + (x^2*asinh(a*x)^(5/2))/2 - (15*sqrt(pi/2)*erf(sqrt(2)*sqrt(asinh(a*x))))/(256*a^2) - (15*sqrt(pi/2)*erfi(sqrt(2)*sqrt(asinh(a*x))))/(256*a^2)]
    @test_int [asinh(a*x)^(5/2), x, 9, (15*x*sqrt(asinh(a*x)))/4 - (5*sqrt(1 + a^2*x^2)*asinh(a*x)^(3/2))/(2*a) + x*asinh(a*x)^(5/2) + (15*sqrt(pi)*erf(sqrt(asinh(a*x))))/(16*a) - (15*sqrt(pi)*erfi(sqrt(asinh(a*x))))/(16*a)]
    @test_int [asinh(a*x)^(5/2)/x, x, 0, Unintegrable(asinh(a*x)^(5/2)/x, x)]

    #= ::Subsubsection::Closed:: =#
    #=n<0=#

    @test_int [x^4/sqrt(asinh(a*x)), x, 18, (sqrt(pi)*erf(sqrt(asinh(a*x))))/(16*a^5) - (sqrt(3*pi)*erf(sqrt(3)*sqrt(asinh(a*x))))/(32*a^5) + (sqrt(pi/5)*erf(sqrt(5)*sqrt(asinh(a*x))))/(32*a^5) + (sqrt(pi)*erfi(sqrt(asinh(a*x))))/(16*a^5) - (sqrt(3*pi)*erfi(sqrt(3)*sqrt(asinh(a*x))))/(32*a^5) + (sqrt(pi/5)*erfi(sqrt(5)*sqrt(asinh(a*x))))/(32*a^5)]
    @test_int [x^3/sqrt(asinh(a*x)), x, 13, -(sqrt(pi)*erf(2*sqrt(asinh(a*x))))/(32*a^4) + (sqrt(pi/2)*erf(sqrt(2)*sqrt(asinh(a*x))))/(8*a^4) + (sqrt(pi)*erfi(2*sqrt(asinh(a*x))))/(32*a^4) - (sqrt(pi/2)*erfi(sqrt(2)*sqrt(asinh(a*x))))/(8*a^4)]
    @test_int [x^2/sqrt(asinh(a*x)), x, 13, -(sqrt(pi)*erf(sqrt(asinh(a*x))))/(8*a^3) + (sqrt(pi/3)*erf(sqrt(3)*sqrt(asinh(a*x))))/(8*a^3) - (sqrt(pi)*erfi(sqrt(asinh(a*x))))/(8*a^3) + (sqrt(pi/3)*erfi(sqrt(3)*sqrt(asinh(a*x))))/(8*a^3)]
    @test_int [x/sqrt(asinh(a*x)), x, 8, -(sqrt(pi/2)*erf(sqrt(2)*sqrt(asinh(a*x))))/(4*a^2) + (sqrt(pi/2)*erfi(sqrt(2)*sqrt(asinh(a*x))))/(4*a^2)]
    @test_int [1/sqrt(asinh(a*x)), x, 6, (sqrt(pi)*erf(sqrt(asinh(a*x))))/(2*a) + (sqrt(pi)*erfi(sqrt(asinh(a*x))))/(2*a)]
    @test_int [1/(x*sqrt(asinh(a*x))), x, 0, Unintegrable(1/(x*sqrt(asinh(a*x))), x)]
    @test_int [1/(x^2*sqrt(asinh(a*x))), x, 0, Unintegrable(1/(x^2*sqrt(asinh(a*x))), x)]

    @test_int [x^4/asinh(a*x)^(3/2), x, 17, -((2*x^4*sqrt(1 + a^2*x^2))/(a*sqrt(asinh(a*x)))) - (sqrt(pi)*erf(sqrt(asinh(a*x))))/(8*a^5) + (3*sqrt(3*pi)*erf(sqrt(3)*sqrt(asinh(a*x))))/(16*a^5) - (sqrt(5*pi)*erf(sqrt(5)*sqrt(asinh(a*x))))/(16*a^5) + (sqrt(pi)*erfi(sqrt(asinh(a*x))))/(8*a^5) - (3*sqrt(3*pi)*erfi(sqrt(3)*sqrt(asinh(a*x))))/(16*a^5) + (sqrt(5*pi)*erfi(sqrt(5)*sqrt(asinh(a*x))))/(16*a^5)]
    @test_int [x^3/asinh(a*x)^(3/2), x, 12, -((2*x^3*sqrt(1 + a^2*x^2))/(a*sqrt(asinh(a*x)))) + (sqrt(pi)*erf(2*sqrt(asinh(a*x))))/(4*a^4) - (sqrt(pi/2)*erf(sqrt(2)*sqrt(asinh(a*x))))/(2*a^4) + (sqrt(pi)*erfi(2*sqrt(asinh(a*x))))/(4*a^4) - (sqrt(pi/2)*erfi(sqrt(2)*sqrt(asinh(a*x))))/(2*a^4)]
    @test_int [x^2/asinh(a*x)^(3/2), x, 12, -((2*x^2*sqrt(1 + a^2*x^2))/(a*sqrt(asinh(a*x)))) + (sqrt(pi)*erf(sqrt(asinh(a*x))))/(4*a^3) - (sqrt(3*pi)*erf(sqrt(3)*sqrt(asinh(a*x))))/(4*a^3) - (sqrt(pi)*erfi(sqrt(asinh(a*x))))/(4*a^3) + (sqrt(3*pi)*erfi(sqrt(3)*sqrt(asinh(a*x))))/(4*a^3)]
    @test_int [x/asinh(a*x)^(3/2), x, 6, (-2*x*sqrt(1 + a^2*x^2))/(a*sqrt(asinh(a*x))) + (sqrt(pi/2)*erf(sqrt(2)*sqrt(asinh(a*x))))/a^2 + (sqrt(pi/2)*erfi(sqrt(2)*sqrt(asinh(a*x))))/a^2]
    @test_int [asinh(a*x)^(-3/2), x, 7, (-2*sqrt(1 + a^2*x^2))/(a*sqrt(asinh(a*x))) - (sqrt(pi)*erf(sqrt(asinh(a*x))))/a + (sqrt(pi)*erfi(sqrt(asinh(a*x))))/a]
    @test_int [1/(x*asinh(a*x)^(3/2)), x, 0, Unintegrable(1/(x*asinh(a*x)^(3/2)), x)]

    @test_int [x^4/asinh(a*x)^(5/2), x, 34, -((2*x^4*sqrt(1 + a^2*x^2))/(3*a*asinh(a*x)^(3/2))) - (16*x^3)/(3*a^2*sqrt(asinh(a*x))) - (20*x^5)/(3*sqrt(asinh(a*x))) + (sqrt(pi)*erf(sqrt(asinh(a*x))))/(12*a^5) - (3*sqrt(3*pi)*erf(sqrt(3)*sqrt(asinh(a*x))))/(8*a^5) + (5*sqrt(5*pi)*erf(sqrt(5)*sqrt(asinh(a*x))))/(24*a^5) + (sqrt(pi)*erfi(sqrt(asinh(a*x))))/(12*a^5) - (3*sqrt(3*pi)*erfi(sqrt(3)*sqrt(asinh(a*x))))/(8*a^5) + (5*sqrt(5*pi)*erfi(sqrt(5)*sqrt(asinh(a*x))))/(24*a^5)]
    @test_int [x^3/asinh(a*x)^(5/2), x, 24, (-2*x^3*sqrt(1 + a^2*x^2))/(3*a*asinh(a*x)^(3/2)) - (4*x^2)/(a^2*sqrt(asinh(a*x))) - (16*x^4)/(3*sqrt(asinh(a*x))) - (2*sqrt(pi)*erf(2*sqrt(asinh(a*x))))/(3*a^4) + (sqrt(2*pi)*erf(sqrt(2)*sqrt(asinh(a*x))))/(3*a^4) + (2*sqrt(pi)*erfi(2*sqrt(asinh(a*x))))/(3*a^4) - (sqrt(2*pi)*erfi(sqrt(2)*sqrt(asinh(a*x))))/(3*a^4)]
    @test_int [x^2/asinh(a*x)^(5/2), x, 22, (-2*x^2*sqrt(1 + a^2*x^2))/(3*a*asinh(a*x)^(3/2)) - (8*x)/(3*a^2*sqrt(asinh(a*x))) - (4*x^3)/sqrt(asinh(a*x)) - (sqrt(pi)*erf(sqrt(asinh(a*x))))/(6*a^3) + (sqrt(3*pi)*erf(sqrt(3)*sqrt(asinh(a*x))))/(2*a^3) - (sqrt(pi)*erfi(sqrt(asinh(a*x))))/(6*a^3) + (sqrt(3*pi)*erfi(sqrt(3)*sqrt(asinh(a*x))))/(2*a^3)]
    @test_int [x/asinh(a*x)^(5/2), x, 11, (-2*x*sqrt(1 + a^2*x^2))/(3*a*asinh(a*x)^(3/2)) - 4/(3*a^2*sqrt(asinh(a*x))) - (8*x^2)/(3*sqrt(asinh(a*x))) - (2*sqrt(2*pi)*erf(sqrt(2)*sqrt(asinh(a*x))))/(3*a^2) + (2*sqrt(2*pi)*erfi(sqrt(2)*sqrt(asinh(a*x))))/(3*a^2)]
    @test_int [1/asinh(a*x)^(5/2), x, 8, (-2*sqrt(1 + a^2*x^2))/(3*a*asinh(a*x)^(3/2)) - (4*x)/(3*sqrt(asinh(a*x))) + (2*sqrt(pi)*erf(sqrt(asinh(a*x))))/(3*a) + (2*sqrt(pi)*erfi(sqrt(asinh(a*x))))/(3*a)]
    @test_int [1/(x*asinh(a*x)^(5/2)), x, 0, Unintegrable(1/(x*asinh(a*x)^(5/2)), x)]

    @test_int [x^4/asinh(a*x)^(7/2), x, 32, -((2*x^4*sqrt(1 + a^2*x^2))/(5*a*asinh(a*x)^(5/2))) - (16*x^3)/(15*a^2*asinh(a*x)^(3/2)) - (4*x^5)/(3*asinh(a*x)^(3/2)) - (32*x^2*sqrt(1 + a^2*x^2))/(5*a^3*sqrt(asinh(a*x))) - (40*x^4*sqrt(1 + a^2*x^2))/(3*a*sqrt(asinh(a*x))) - (sqrt(pi)*erf(sqrt(asinh(a*x))))/(30*a^5) + (9*sqrt(3*pi)*erf(sqrt(3)*sqrt(asinh(a*x))))/(20*a^5) - (5*sqrt(5*pi)*erf(sqrt(5)*sqrt(asinh(a*x))))/(12*a^5) + (sqrt(pi)*erfi(sqrt(asinh(a*x))))/(30*a^5) - (9*sqrt(3*pi)*erfi(sqrt(3)*sqrt(asinh(a*x))))/(20*a^5) + (5*sqrt(5*pi)*erfi(sqrt(5)*sqrt(asinh(a*x))))/(12*a^5)]
    @test_int [x^3/asinh(a*x)^(7/2), x, 21, (-2*x^3*sqrt(1 + a^2*x^2))/(5*a*asinh(a*x)^(5/2)) - (4*x^2)/(5*a^2*asinh(a*x)^(3/2)) - (16*x^4)/(15*asinh(a*x)^(3/2)) - (16*x*sqrt(1 + a^2*x^2))/(5*a^3*sqrt(asinh(a*x))) - (128*x^3*sqrt(1 + a^2*x^2))/(15*a*sqrt(asinh(a*x))) + (16*sqrt(pi)*erf(2*sqrt(asinh(a*x))))/(15*a^4) - (4*sqrt(2*pi)*erf(sqrt(2)*sqrt(asinh(a*x))))/(15*a^4) + (16*sqrt(pi)*erfi(2*sqrt(asinh(a*x))))/(15*a^4) - (4*sqrt(2*pi)*erfi(sqrt(2)*sqrt(asinh(a*x))))/(15*a^4)]
    @test_int [x^2/asinh(a*x)^(7/2), x, 22, (-2*x^2*sqrt(1 + a^2*x^2))/(5*a*asinh(a*x)^(5/2)) - (8*x)/(15*a^2*asinh(a*x)^(3/2)) - (4*x^3)/(5*asinh(a*x)^(3/2)) - (16*sqrt(1 + a^2*x^2))/(15*a^3*sqrt(asinh(a*x))) - (24*x^2*sqrt(1 + a^2*x^2))/(5*a*sqrt(asinh(a*x))) + (sqrt(pi)*erf(sqrt(asinh(a*x))))/(15*a^3) - (3*sqrt(3*pi)*erf(sqrt(3)*sqrt(asinh(a*x))))/(5*a^3) - (sqrt(pi)*erfi(sqrt(asinh(a*x))))/(15*a^3) + (3*sqrt(3*pi)*erfi(sqrt(3)*sqrt(asinh(a*x))))/(5*a^3)]
    @test_int [x/asinh(a*x)^(7/2), x, 9, (-2*x*sqrt(1 + a^2*x^2))/(5*a*asinh(a*x)^(5/2)) - 4/(15*a^2*asinh(a*x)^(3/2)) - (8*x^2)/(15*asinh(a*x)^(3/2)) - (32*x*sqrt(1 + a^2*x^2))/(15*a*sqrt(asinh(a*x))) + (8*sqrt(2*pi)*erf(sqrt(2)*sqrt(asinh(a*x))))/(15*a^2) + (8*sqrt(2*pi)*erfi(sqrt(2)*sqrt(asinh(a*x))))/(15*a^2)]
    @test_int [asinh(a*x)^(-7/2), x, 9, (-2*sqrt(1 + a^2*x^2))/(5*a*asinh(a*x)^(5/2)) - (4*x)/(15*asinh(a*x)^(3/2)) - (8*sqrt(1 + a^2*x^2))/(15*a*sqrt(asinh(a*x))) - (4*sqrt(pi)*erf(sqrt(asinh(a*x))))/(15*a) + (4*sqrt(pi)*erfi(sqrt(asinh(a*x))))/(15*a)]
    @test_int [1/(x*asinh(a*x)^(7/2)), x, 0, Unintegrable(1/(x*asinh(a*x)^(7/2)), x)]

    #= ::Subsection::Closed:: =#
    #=Integrands*of*the*form*(b*x)^m*asin(a*x)^n*with*m*symbolic=#

    @test_int [x^m*asinh(a*x)^4, x, 1, (x^(1 + m)*asinh(a*x)^4)/(1 + m) - (4*a*Unintegrable((x^(1 + m)*asinh(a*x)^3)/sqrt(1 + a^2*x^2), x))/(1 + m)]
    @test_int [x^m*asinh(a*x)^3, x, 1, (x^(1 + m)*asinh(a*x)^3)/(1 + m) - (3*a*Unintegrable((x^(1 + m)*asinh(a*x)^2)/sqrt(1 + a^2*x^2), x))/(1 + m)]
    @test_int [x^m*asinh(a*x)^2, x, 2, (x^(1 + m)*asinh(a*x)^2)/(1 + m) - (2*a*x^(2 + m)*asinh(a*x)*HypergeometricFunctions._₂F₁(1/2, (2 + m)/2, (4 + m)/2, (-a^2)*x^2))/(2 + 3*m + m^2) + (2*a^2*x^(3 + m)*HypergeometricPFQ([1, 3/2 + m/2, 3/2 + m/2], [2 + m/2, 5/2 + m/2], (-a^2)*x^2))/(6 + 11*m + 6*m^2 + m^3)]
    @test_int [x^m*asinh(a*x)^1, x, 2, (x^(1 + m)*asinh(a*x))/(1 + m) - (a*x^(2 + m)*HypergeometricFunctions._₂F₁(1/2, (2 + m)/2, (4 + m)/2, -(a^2*x^2)))/(2 + 3*m + m^2)]
    @test_int [x^m/asinh(a*x)^1, x, 0, Unintegrable(x^m/asinh(a*x), x)]
    @test_int [x^m/asinh(a*x)^2, x, 0, Unintegrable(x^m/asinh(a*x)^2, x)]

    @test_int [x^m*asinh(a*x)^(5/2), x, 0, Unintegrable(x^m*asinh(a*x)^(5/2), x)]
    @test_int [x^m*asinh(a*x)^(3/2), x, 0, Unintegrable(x^m*asinh(a*x)^(3/2), x)]
    @test_int [x^m*sqrt(asinh(a*x)), x, 0, Unintegrable(x^m*sqrt(asinh(a*x)), x)]
    @test_int [x^m/sqrt(asinh(a*x)), x, 0, Unintegrable(x^m/sqrt(asinh(a*x)), x)]
    @test_int [x^m/asinh(a*x)^(3/2), x, 0, Unintegrable(x^m/asinh(a*x)^(3/2), x)]

    #= ::Subsection::Closed:: =#
    #=Integrands*of*the*form*(b*x)^m*asin(a*x)^n*with*n*symbolic=#

    @test_int [(b*x)^m*asinh(a*x)^n, x, 0, Unintegrable((b*x)^m*asinh(a*x)^n, x)]

    @test_int [x^4*asinh(a*x)^n, x, 12, (5^(-1 - n)*asinh(a*x)^n*Gamma(1 + n, -5*asinh(a*x)))/((-asinh(a*x))^n*(32*a^5)) - (asinh(a*x)^n*Gamma(1 + n, -3*asinh(a*x)))/(3^n*(-asinh(a*x))^n*(32*a^5)) + (asinh(a*x)^n*Gamma(1 + n, -asinh(a*x)))/((-asinh(a*x))^n*(16*a^5)) - Gamma(1 + n, asinh(a*x))/(16*a^5) + Gamma(1 + n, 3*asinh(a*x))/(3^n*(32*a^5)) - (5^(-1 - n)*Gamma(1 + n, 5*asinh(a*x)))/(32*a^5)]
    @test_int [x^3*asinh(a*x)^n, x, 9, (asinh(a*x)^n*Gamma(1 + n, -4*asinh(a*x)))/(2^(2*(3 + n))*(-asinh(a*x))^n*a^4) - (2^(-4 - n)*asinh(a*x)^n*Gamma(1 + n, -2*asinh(a*x)))/((-asinh(a*x))^n*a^4) - (2^(-4 - n)*Gamma(1 + n, 2*asinh(a*x)))/a^4 + Gamma(1 + n, 4*asinh(a*x))/(2^(2*(3 + n))*a^4)]
    @test_int [x^2*asinh(a*x)^n, x, 9, (3^(-1 - n)*asinh(a*x)^n*Gamma(1 + n, -3*asinh(a*x)))/((-asinh(a*x))^n*(8*a^3)) - (asinh(a*x)^n*Gamma(1 + n, -asinh(a*x)))/((-asinh(a*x))^n*(8*a^3)) + Gamma(1 + n, asinh(a*x))/(8*a^3) - (3^(-1 - n)*Gamma(1 + n, 3*asinh(a*x)))/(8*a^3)]
    @test_int [x^1*asinh(a*x)^n, x, 6, (2^(-3 - n)*asinh(a*x)^n*Gamma(1 + n, -2*asinh(a*x)))/((-asinh(a*x))^n*a^2) + (2^(-3 - n)*Gamma(1 + n, 2*asinh(a*x)))/a^2]
    @test_int [x^0*asinh(a*x)^n, x, 4, (asinh(a*x)^n*Gamma(1 + n, -asinh(a*x)))/(2*a*(-asinh(a*x))^n) - Gamma(1 + n, asinh(a*x))/(2*a)]
    @test_int [asinh(a*x)^n/x^1, x, 0, Unintegrable(asinh(a*x)^n/x, x)]
    @test_int [asinh(a*x)^n/x^2, x, 0, Unintegrable(asinh(a*x)^n/x^2, x)]

    #= ::Section::Closed:: =#
    #=Integrands*of*the*form*(d*x)^m*(a+b*asinh(c*x)^p=#

    #= ::Subsection::Closed:: =#
    #=Integrands*of*the*form*x^m*(a+b*asinh(c*x))^(n/2)=#

    #= ::Subsubsection::Closed:: =#
    #=n>0=#

    @test_int [x^2*(a + b*asinh(c*x))^(1/2), x, 14, (1/3)*x^3*sqrt(a + b*asinh(c*x)) - (sqrt(b)*ℯ^(a/b)*sqrt(pi)*erf(sqrt(a + b*asinh(c*x))/sqrt(b)))/(16*c^3) + (sqrt(b)*ℯ^((3*a)/b)*sqrt(pi/3)*erf((sqrt(3)*sqrt(a + b*asinh(c*x)))/sqrt(b)))/(48*c^3) + (sqrt(b)*sqrt(pi)*erfi(sqrt(a + b*asinh(c*x))/sqrt(b)))/(ℯ^(a/b)*(16*c^3)) - (sqrt(b)*sqrt(pi/3)*erfi((sqrt(3)*sqrt(a + b*asinh(c*x)))/sqrt(b)))/(ℯ^((3*a)/b)*(48*c^3))]
    @test_int [x^1*(a + b*asinh(c*x))^(1/2), x, 9, sqrt(a + b*asinh(c*x))/(4*c^2) + (1/2)*x^2*sqrt(a + b*asinh(c*x)) - (sqrt(b)*ℯ^((2*a)/b)*sqrt(pi/2)*erf((sqrt(2)*sqrt(a + b*asinh(c*x)))/sqrt(b)))/(16*c^2) - (sqrt(b)*sqrt(pi/2)*erfi((sqrt(2)*sqrt(a + b*asinh(c*x)))/sqrt(b)))/(ℯ^((2*a)/b)*(16*c^2))]
    @test_int [x^0*(a + b*asinh(c*x))^(1/2), x, 7, x*sqrt(a + b*asinh(c*x)) + (sqrt(b)*ℯ^(a/b)*sqrt(pi)*erf(sqrt(a + b*asinh(c*x))/sqrt(b)))/(4*c) - (sqrt(b)*sqrt(pi)*erfi(sqrt(a + b*asinh(c*x))/sqrt(b)))/(ℯ^(a/b)*(4*c))]

    @test_int [x^2*(a + b*asinh(c*x))^(3/2), x, 22, (b*sqrt(1 + c^2*x^2)*sqrt(a + b*asinh(c*x)))/(3*c^3) - (b*x^2*sqrt(1 + c^2*x^2)*sqrt(a + b*asinh(c*x)))/(6*c) + (1/3)*x^3*(a + b*asinh(c*x))^(3/2) - (3*b^(3/2)*ℯ^(a/b)*sqrt(pi)*erf(sqrt(a + b*asinh(c*x))/sqrt(b)))/(32*c^3) + (b^(3/2)*ℯ^((3*a)/b)*sqrt(pi/3)*erf((sqrt(3)*sqrt(a + b*asinh(c*x)))/sqrt(b)))/(96*c^3) - (3*b^(3/2)*sqrt(pi)*erfi(sqrt(a + b*asinh(c*x))/sqrt(b)))/(ℯ^(a/b)*(32*c^3)) + (b^(3/2)*sqrt(pi/3)*erfi((sqrt(3)*sqrt(a + b*asinh(c*x)))/sqrt(b)))/(ℯ^((3*a)/b)*(96*c^3))]
    @test_int [x^1*(a + b*asinh(c*x))^(3/2), x, 11, -((3*b*x*sqrt(1 + c^2*x^2)*sqrt(a + b*asinh(c*x)))/(8*c)) + (a + b*asinh(c*x))^(3/2)/(4*c^2) + (1/2)*x^2*(a + b*asinh(c*x))^(3/2) - (3*b^(3/2)*ℯ^((2*a)/b)*sqrt(pi/2)*erf((sqrt(2)*sqrt(a + b*asinh(c*x)))/sqrt(b)))/(64*c^2) + (3*b^(3/2)*sqrt(pi/2)*erfi((sqrt(2)*sqrt(a + b*asinh(c*x)))/sqrt(b)))/(ℯ^((2*a)/b)*(64*c^2))]
    @test_int [x^0*(a + b*asinh(c*x))^(3/2), x, 8, -((3*b*sqrt(1 + c^2*x^2)*sqrt(a + b*asinh(c*x)))/(2*c)) + x*(a + b*asinh(c*x))^(3/2) + (3*b^(3/2)*ℯ^(a/b)*sqrt(pi)*erf(sqrt(a + b*asinh(c*x))/sqrt(b)))/(8*c) + (3*b^(3/2)*sqrt(pi)*erfi(sqrt(a + b*asinh(c*x))/sqrt(b)))/(ℯ^(a/b)*(8*c))]

    @test_int [x^2*(a + b*asinh(c*x))^(5/2), x, 24, -((5*b^2*x*sqrt(a + b*asinh(c*x)))/(6*c^2)) + (5/36)*b^2*x^3*sqrt(a + b*asinh(c*x)) + (5*b*sqrt(1 + c^2*x^2)*(a + b*asinh(c*x))^(3/2))/(9*c^3) - (5*b*x^2*sqrt(1 + c^2*x^2)*(a + b*asinh(c*x))^(3/2))/(18*c) + (1/3)*x^3*(a + b*asinh(c*x))^(5/2) - (15*b^(5/2)*ℯ^(a/b)*sqrt(pi)*erf(sqrt(a + b*asinh(c*x))/sqrt(b)))/(64*c^3) + (5*b^(5/2)*ℯ^((3*a)/b)*sqrt(pi/3)*erf((sqrt(3)*sqrt(a + b*asinh(c*x)))/sqrt(b)))/(576*c^3) + (15*b^(5/2)*sqrt(pi)*erfi(sqrt(a + b*asinh(c*x))/sqrt(b)))/(ℯ^(a/b)*(64*c^3)) - (5*b^(5/2)*sqrt(pi/3)*erfi((sqrt(3)*sqrt(a + b*asinh(c*x)))/sqrt(b)))/(ℯ^((3*a)/b)*(576*c^3))]
    @test_int [x^1*(a + b*asinh(c*x))^(5/2), x, 12, (15*b^2*sqrt(a + b*asinh(c*x)))/(64*c^2) + (15/32)*b^2*x^2*sqrt(a + b*asinh(c*x)) - (5*b*x*sqrt(1 + c^2*x^2)*(a + b*asinh(c*x))^(3/2))/(8*c) + (a + b*asinh(c*x))^(5/2)/(4*c^2) + (1/2)*x^2*(a + b*asinh(c*x))^(5/2) - (15*b^(5/2)*ℯ^((2*a)/b)*sqrt(pi/2)*erf((sqrt(2)*sqrt(a + b*asinh(c*x)))/sqrt(b)))/(256*c^2) - (15*b^(5/2)*sqrt(pi/2)*erfi((sqrt(2)*sqrt(a + b*asinh(c*x)))/sqrt(b)))/(ℯ^((2*a)/b)*(256*c^2))]
    @test_int [x^0*(a + b*asinh(c*x))^(5/2), x, 9, (15/4)*b^2*x*sqrt(a + b*asinh(c*x)) - (5*b*sqrt(1 + c^2*x^2)*(a + b*asinh(c*x))^(3/2))/(2*c) + x*(a + b*asinh(c*x))^(5/2) + (15*b^(5/2)*ℯ^(a/b)*sqrt(pi)*erf(sqrt(a + b*asinh(c*x))/sqrt(b)))/(16*c) - (15*b^(5/2)*sqrt(pi)*erfi(sqrt(a + b*asinh(c*x))/sqrt(b)))/(ℯ^(a/b)*(16*c))]

    #= ::Subsubsection::Closed:: =#
    #=n<0=#

    @test_int [x^2/(a + b*asinh(c*x))^(1/2), x, 13, -((ℯ^(a/b)*sqrt(pi)*erf(sqrt(a + b*asinh(c*x))/sqrt(b)))/(8*sqrt(b)*c^3)) + (ℯ^((3*a)/b)*sqrt(pi/3)*erf((sqrt(3)*sqrt(a + b*asinh(c*x)))/sqrt(b)))/(8*sqrt(b)*c^3) - (sqrt(pi)*erfi(sqrt(a + b*asinh(c*x))/sqrt(b)))/(ℯ^(a/b)*(8*sqrt(b)*c^3)) + (sqrt(pi/3)*erfi((sqrt(3)*sqrt(a + b*asinh(c*x)))/sqrt(b)))/(ℯ^((3*a)/b)*(8*sqrt(b)*c^3))]
    @test_int [x^1/(a + b*asinh(c*x))^(1/2), x, 8, -((ℯ^((2*a)/b)*sqrt(pi/2)*erf((sqrt(2)*sqrt(a + b*asinh(c*x)))/sqrt(b)))/(4*sqrt(b)*c^2)) + (sqrt(pi/2)*erfi((sqrt(2)*sqrt(a + b*asinh(c*x)))/sqrt(b)))/(ℯ^((2*a)/b)*(4*sqrt(b)*c^2))]
    @test_int [x^0/(a + b*asinh(c*x))^(1/2), x, 6, (ℯ^(a/b)*sqrt(pi)*erf(sqrt(a + b*asinh(c*x))/sqrt(b)))/(2*sqrt(b)*c) + (sqrt(pi)*erfi(sqrt(a + b*asinh(c*x))/sqrt(b)))/(ℯ^(a/b)*(2*sqrt(b)*c))]

    @test_int [x^2/(a + b*asinh(c*x))^(3/2), x, 12, -((2*x^2*sqrt(1 + c^2*x^2))/(b*c*sqrt(a + b*asinh(c*x)))) + (ℯ^(a/b)*sqrt(pi)*erf(sqrt(a + b*asinh(c*x))/sqrt(b)))/(4*b^(3/2)*c^3) - (ℯ^((3*a)/b)*sqrt(3*pi)*erf((sqrt(3)*sqrt(a + b*asinh(c*x)))/sqrt(b)))/(4*b^(3/2)*c^3) - (sqrt(pi)*erfi(sqrt(a + b*asinh(c*x))/sqrt(b)))/(ℯ^(a/b)*(4*b^(3/2)*c^3)) + (sqrt(3*pi)*erfi((sqrt(3)*sqrt(a + b*asinh(c*x)))/sqrt(b)))/(ℯ^((3*a)/b)*(4*b^(3/2)*c^3))]
    @test_int [x^1/(a + b*asinh(c*x))^(3/2), x, 6, -((2*x*sqrt(1 + c^2*x^2))/(b*c*sqrt(a + b*asinh(c*x)))) + (ℯ^((2*a)/b)*sqrt(pi/2)*erf((sqrt(2)*sqrt(a + b*asinh(c*x)))/sqrt(b)))/(b^(3/2)*c^2) + (sqrt(pi/2)*erfi((sqrt(2)*sqrt(a + b*asinh(c*x)))/sqrt(b)))/(ℯ^((2*a)/b)*(b^(3/2)*c^2))]
    @test_int [x^0/(a + b*asinh(c*x))^(3/2), x, 7, -((2*sqrt(1 + c^2*x^2))/(b*c*sqrt(a + b*asinh(c*x)))) - (ℯ^(a/b)*sqrt(pi)*erf(sqrt(a + b*asinh(c*x))/sqrt(b)))/(b^(3/2)*c) + (sqrt(pi)*erfi(sqrt(a + b*asinh(c*x))/sqrt(b)))/(ℯ^(a/b)*(b^(3/2)*c))]

    @test_int [x^2/(a + b*asinh(c*x))^(5/2), x, 22, -((2*x^2*sqrt(1 + c^2*x^2))/(3*b*c*(a + b*asinh(c*x))^(3/2))) - (8*x)/(3*b^2*c^2*sqrt(a + b*asinh(c*x))) - (4*x^3)/(b^2*sqrt(a + b*asinh(c*x))) - (ℯ^(a/b)*sqrt(pi)*erf(sqrt(a + b*asinh(c*x))/sqrt(b)))/(6*b^(5/2)*c^3) + (ℯ^((3*a)/b)*sqrt(3*pi)*erf((sqrt(3)*sqrt(a + b*asinh(c*x)))/sqrt(b)))/(2*b^(5/2)*c^3) - (sqrt(pi)*erfi(sqrt(a + b*asinh(c*x))/sqrt(b)))/(ℯ^(a/b)*(6*b^(5/2)*c^3)) + (sqrt(3*pi)*erfi((sqrt(3)*sqrt(a + b*asinh(c*x)))/sqrt(b)))/(ℯ^((3*a)/b)*(2*b^(5/2)*c^3))]
    @test_int [x^1/(a + b*asinh(c*x))^(5/2), x, 11, -((2*x*sqrt(1 + c^2*x^2))/(3*b*c*(a + b*asinh(c*x))^(3/2))) - 4/(3*b^2*c^2*sqrt(a + b*asinh(c*x))) - (8*x^2)/(3*b^2*sqrt(a + b*asinh(c*x))) - (2*ℯ^((2*a)/b)*sqrt(2*pi)*erf((sqrt(2)*sqrt(a + b*asinh(c*x)))/sqrt(b)))/(3*b^(5/2)*c^2) + (2*sqrt(2*pi)*erfi((sqrt(2)*sqrt(a + b*asinh(c*x)))/sqrt(b)))/(ℯ^((2*a)/b)*(3*b^(5/2)*c^2))]
    @test_int [x^0/(a + b*asinh(c*x))^(5/2), x, 8, -((2*sqrt(1 + c^2*x^2))/(3*b*c*(a + b*asinh(c*x))^(3/2))) - (4*x)/(3*b^2*sqrt(a + b*asinh(c*x))) + (2*ℯ^(a/b)*sqrt(pi)*erf(sqrt(a + b*asinh(c*x))/sqrt(b)))/(3*b^(5/2)*c) + (2*sqrt(pi)*erfi(sqrt(a + b*asinh(c*x))/sqrt(b)))/(ℯ^(a/b)*(3*b^(5/2)*c))]

    @test_int [x^2/(a + b*asinh(c*x))^(7/2), x, 22, -((2*x^2*sqrt(1 + c^2*x^2))/(5*b*c*(a + b*asinh(c*x))^(5/2))) - (8*x)/(15*b^2*c^2*(a + b*asinh(c*x))^(3/2)) - (4*x^3)/(5*b^2*(a + b*asinh(c*x))^(3/2)) - (16*sqrt(1 + c^2*x^2))/(15*b^3*c^3*sqrt(a + b*asinh(c*x))) - (24*x^2*sqrt(1 + c^2*x^2))/(5*b^3*c*sqrt(a + b*asinh(c*x))) + (ℯ^(a/b)*sqrt(pi)*erf(sqrt(a + b*asinh(c*x))/sqrt(b)))/(15*b^(7/2)*c^3) - (3*ℯ^((3*a)/b)*sqrt(3*pi)*erf((sqrt(3)*sqrt(a + b*asinh(c*x)))/sqrt(b)))/(5*b^(7/2)*c^3) - (sqrt(pi)*erfi(sqrt(a + b*asinh(c*x))/sqrt(b)))/(ℯ^(a/b)*(15*b^(7/2)*c^3)) + (3*sqrt(3*pi)*erfi((sqrt(3)*sqrt(a + b*asinh(c*x)))/sqrt(b)))/(ℯ^((3*a)/b)*(5*b^(7/2)*c^3))]
    @test_int [x^1/(a + b*asinh(c*x))^(7/2), x, 9, -((2*x*sqrt(1 + c^2*x^2))/(5*b*c*(a + b*asinh(c*x))^(5/2))) - 4/(15*b^2*c^2*(a + b*asinh(c*x))^(3/2)) - (8*x^2)/(15*b^2*(a + b*asinh(c*x))^(3/2)) - (32*x*sqrt(1 + c^2*x^2))/(15*b^3*c*sqrt(a + b*asinh(c*x))) + (8*ℯ^((2*a)/b)*sqrt(2*pi)*erf((sqrt(2)*sqrt(a + b*asinh(c*x)))/sqrt(b)))/(15*b^(7/2)*c^2) + (8*sqrt(2*pi)*erfi((sqrt(2)*sqrt(a + b*asinh(c*x)))/sqrt(b)))/(ℯ^((2*a)/b)*(15*b^(7/2)*c^2))]
    @test_int [x^0/(a + b*asinh(c*x))^(7/2), x, 9, -((2*sqrt(1 + c^2*x^2))/(5*b*c*(a + b*asinh(c*x))^(5/2))) - (4*x)/(15*b^2*(a + b*asinh(c*x))^(3/2)) - (8*sqrt(1 + c^2*x^2))/(15*b^3*c*sqrt(a + b*asinh(c*x))) - (4*ℯ^(a/b)*sqrt(pi)*erf(sqrt(a + b*asinh(c*x))/sqrt(b)))/(15*b^(7/2)*c) + (4*sqrt(pi)*erfi(sqrt(a + b*asinh(c*x))/sqrt(b)))/(ℯ^(a/b)*(15*b^(7/2)*c))]
    end

max_stars_repo_path: test/7 Inverse hyperbolic functions/7.1 Inverse hyperbolic sine/7.1.2 (d x)^m (a+b arcsinh(c x))^n.jl


text:
# Used internally to expose open disk files inside a `do` block
struct OpenGeoArray{T,N,D<:Tuple,R<:Tuple,A,Na<:Symbol,Me,Mi} <: AbstractGeoArray{T,N,D,A}
    data::A
    dims::D
    refdims::R
    name::Na
    metadata::Me
    missingval::Mi
end
function OpenGeoArray(A, dims, refdims, name, metadata, missingval)
    OpenGeoArray{eltype(A),ndims(A),map(typeof,(dims,refdims,A,name,metadata,missingval))...}(
        A, dims, refdims, name, metadata, missingval
    )
end
function OpenGeoArray(f, A::AbstractGeoArray{T,N}) where {T,N}
    withsourcedata(A) do source
        OA = OpenGeoArray(source, dims(A), refdims(A), name(A), metadata(A), missingval(A))
        f(OA)
    end
end

"""
    Open(f, A::AbstractGeoArray)

`Open` is used to open any `AbstractGeoArray` and do multiple operations
on it in a safe way. It's a shorthand for the unexported `OpenGeoArray`
constructor.

`f` is a method that accepts a single argument - an `OpenGeoArray` object
which is just an `AbstractGeoArray` that holds an open disk - based object.
Often it will be a `do` block:

```julia
ga = GDALarray(filepath)
Open(ga) do A
    A[I...] # A is an `OpenGeoArray` wrapping the disk-based object.
    # ...  multiple things you need to do with the open file
end
```

By using a do block to open file we ensure they are always closed again
after we finish working with them.
"""
Open(f, A::AbstractGeoArray) = OpenGeoArray(f, A)


text:
# -*- coding: utf-8 -*-
# ---
# jupyter:
#   jupytext:
#     formats: ipynb,jl:light
#     text_representation:
#       extension: .jl
#       format_name: light
#       format_version: '1.4'
#       jupytext_version: 1.2.4
#   kernelspec:
#     display_name: Julia 1.3.1
#     language: julia
#     name: julia-1.3
# ---

# + {"slideshow": {"slide_type": "slide"}, "cell_type": "markdown"}
# # Lagrangian particle tracking
#
# Material particles following ocean currents can be analyzed in terms of trajectories. These can simply be computed by integrating velocities through time within a [Lagrangian framework](https://en.wikipedia.org/wiki/Lagrangian_and_Eulerian_specification_of_the_flow_field). In `Julia` this is easily done using the [IndividualDisplacements.jl](https://JuliaClimate.github.io/IndividualDisplacements.jl/dev/) and [OrdinaryDiffEq.jl](https://docs.juliadiffeq.org/latest) packages.

# + {"slideshow": {"slide_type": "slide"}, "cell_type": "markdown"}
# ## Import software, grid, and velocity fields
#
# 1. pre-requisites
# 2. read variables

# + {"slideshow": {"slide_type": "subslide"}}
#]add MITgcmTools#master; add OrdinaryDiffEq; add IndividualDisplacements#master;

# + {"slideshow": {"slide_type": "-"}, "cell_style": "center"}
using IndividualDisplacements, MeshArrays, OrdinaryDiffEq
using Plots, Statistics, MITgcmTools, DataFrames

include("helper_functions.jl")
get_grid_if_needed()
get_velocity_if_needed()

Γ=read_llc90_grid();

# + {"slideshow": {"slide_type": "slide"}, "cell_type": "markdown"}
# ## Set Gridded Variables
#
# 1. Select depth
# 2. time average & normalize to grid units
# 3. apply exchange function to `u,v,lon,lat`
# 4. store everything in `uv_etc` dictionary

# + {"slideshow": {"slide_type": "-"}}
uv_etc=read_uv_etc(20,Γ);

# + {"slideshow": {"slide_type": "slide"}, "cell_type": "markdown"}
# ## Initialize Trajectory Computation
#
# 1. set `duComp` to interpolation function
# 2. set `u0` initial location array

# + {"slideshow": {"slide_type": "-"}}
du_dt=IndividualDisplacements.VelComp!
(u0,du)=initialize_locations(uv_etc,10);

# + {"slideshow": {"slide_type": "slide"}, "cell_type": "markdown"}
# ## Compute Trajectories
#
# _Note: `ODEProblem` and `solve` settings can still be refined_
# -

tspan = (0.0,uv_etc["t1"]-uv_etc["t0"])
prob = ODEProblem(du_dt,u0,tspan,uv_etc)
sol = solve(prob,Euler(),dt=uv_etc["dt"])
size(sol)

# + {"slideshow": {"slide_type": "slide"}, "cell_type": "markdown"}
# ## Post-Process Trajectories
#
# 1. Copy `sol` to a `DataFrame`
# 2. Map position to lon,lat coordinates

# + {"slideshow": {"slide_type": "subslide"}}
df=postprocess_ODESolution(sol,uv_etc);

# + {"slideshow": {"slide_type": "slide"}, "cell_type": "markdown"}
# ## Plot Trajectories
#
#
# For example, you can:
#
# - use either `Makie.jl` via `PlotMakie`
# - or use `PyPlot.jl` via `PlotMapProj`
#

# + {"slideshow": {"slide_type": "-"}}
p=dirname(pathof(IndividualDisplacements));

# + {"slideshow": {"slide_type": "slide"}}
#include(joinpath(p,"plot_pyplot.jl"));
#PyPlot.figure(); PlotMapProj(df,5000)

#include(joinpath(p,"plot_makie.jl")); AbstractPlotting.inline!(true); #for Juno, set to false
#scene=PlotMakie(df,5000,180.0) #Makie.save("LatLonCap300mDepth.png", scene)

include(joinpath(p,"plot_plots.jl"));
plt=PlotBasic(df,1000,180.0)

max_stars_repo_name: abelolivelanna/GlobalOceanNotebooks

max_issues_repo_name: abelolivelanna/GlobalOceanNotebooks

max_forks_repo_name: abelolivelanna/GlobalOceanNotebooks


text:
"""
    FullSearch(d=Cityblock(), N; patch_radius, search_radius, search_stride=1)

Construct a block matching algorithm with full search strategy. In some literature, this is also
called _exhaustive search_.

The default patch distance is `Distances.Cityblock()`.

!!! note
    If `p` is not given and if inputs are `CuArray` then GPU methods will be applied for a small
    set of predefined distances from Distances.jl. The result will be slightly different from its
    CPU version.
"""
struct FullSearch{F,N} <: AbstractBlockMatchingStrategy
    f::F
    patch_radius::CartesianIndex{N}
    search_radius::CartesianIndex{N}
    search_stride::CartesianIndex{N}
end

FullSearch(N; kwargs...) = FullSearch(Cityblock(), N; kwargs...)
FullSearch(f, N; patch_radius, search_radius, search_stride=1) =
    FullSearch{typeof(f),N}(f, to_cartesian(N, patch_radius, search_radius, search_stride)...)

function best_match(S::FullSearch, ref, frame, p::CartesianIndex; offset=false)
    size_check(ref, frame)

    rₚ = S.patch_radius
    R_frame = CartesianIndices(frame)
    R_frame = first(R_frame) + rₚ : last(R_frame) - rₚ
    p in R_frame || throw(ArgumentError("Boundary unsupported now: pixel `p = $p` should be within $(first(R_frame)):$(last(R_frame))."))

    R_ref = CartesianIndices(ref)
    R_ref = first(R_ref) + rₚ : last(R_ref) - rₚ
    q_start = max(first(R_ref), p - S.search_radius)
    q_stop = min(last(R_ref), p + S.search_radius)
    candidates = q_start:S.search_stride:q_stop

    patch_p = frame[p - S.patch_radius:p + S.patch_radius]
    dist = map(candidates) do q
        patch_q = @view ref[q - S.patch_radius:q + S.patch_radius]
        S.f(patch_p, patch_q)
    end

    _, idx = findmin(dist)
    q = candidates[idx]

    return offset ? q - p : q
end

function best_match(S::FullSearch, ref, frame; offset=false)
    size_check(ref, frame)

    rₚ = S.patch_radius

    R_frame = CartesianIndices(frame)
    R_ref = CartesianIndices(ref)
    # does not consider boundary
    R_frame = first(R_frame) + rₚ : last(R_frame) - rₚ
    R_ref = first(R_ref) + rₚ : last(R_ref) - rₚ
    matches = OffsetArray(similar(CartesianIndices(R_frame)), rₚ.I)

    # pre-allocation to reduce memeory allocation
    p = CartesianIndex((last(R_frame) + first(R_frame)).I .÷ 2) # use center point to initialize
    q_start = max(first(R_ref), p - S.search_radius)
    q_stop = min(last(R_ref), p + S.search_radius)
    candidates = q_start:S.search_stride:q_stop
    @assert length(candidates) > 0
    patch_p = frame[p-rₚ:p+rₚ] # pre-allocation into a contiguous memory layout
    dist = Vector{Float32}(undef, length(candidates))

    for p in R_frame
        patch_p .= @view frame[p-rₚ:p+rₚ]

        q_start = max(first(R_ref), p - S.search_radius)
        q_stop = min(last(R_ref), p + S.search_radius)
        candidates = q_start:S.search_stride:q_stop
        n = length(candidates)
        for k in 1:n
            q = candidates[k]
            patch_q = @view ref[q-rₚ:q+rₚ]
            dist[k] = S.f(patch_p, patch_q)
        end
        _, idx = findmin(view(dist, 1:n))
        matches[p] = candidates[idx]
    end

    offset && (matches .= matches .- OffsetArray(R_frame, rₚ.I))
    return matches
end

function multi_match(S::FullSearch, ref, frame, p::CartesianIndex; num_patches, offset=false)
    size_check(ref, frame)

    rₚ = S.patch_radius
    R_frame = CartesianIndices(frame)
    R_frame = first(R_frame) + rₚ : last(R_frame) - rₚ
    p in R_frame || throw(ArgumentError("Boundary unsupported now: pixel `p = $p` should be within $(first(R_frame)):$(last(R_frame))."))

    R_ref = CartesianIndices(ref)
    R_ref = first(R_ref) + rₚ : last(R_ref) - rₚ
    q_start = max(first(R_ref), p - S.search_radius)
    q_stop = min(last(R_ref), p + S.search_radius)
    candidates = q_start:S.search_stride:q_stop

    patch_p = frame[p - S.patch_radius:p + S.patch_radius]
    T = typeof(S.f(patch_p, patch_p))
    dist = Vector{T}(undef, length(candidates))
    @inbounds for i in 1:length(dist)
        q = candidates[i]
        patch_q = @view ref[q - S.patch_radius:q + S.patch_radius]
        dist[i] = S.f(patch_p, patch_q)
    end

    qs = candidates[partialsortperm(dist, 1:num_patches)]
    return offset ? map(q->q-p, qs) : qs
end

function multi_match(S::FullSearch, ref, frame; num_patches, offset=false)
    size_check(ref, frame)

    rₚ = S.patch_radius

    R_frame = CartesianIndices(frame)
    R_ref = CartesianIndices(ref)
    # does not consider boundary
    R_frame = first(R_frame) + rₚ : last(R_frame) - rₚ
    R_ref = first(R_ref) + rₚ : last(R_ref) - rₚ

    matches = Array{CartesianIndex{ndims(ref)}}(undef, size(R_frame)..., num_patches)

    # pre-allocation to reduce memeory allocation
    p = CartesianIndex((last(R_frame) + first(R_frame)).I .÷ 2) # use center point to initialize
    q_start = max(first(R_ref), p - S.search_radius)
    q_stop = min(last(R_ref), p + S.search_radius)
    candidates = q_start:S.search_stride:q_stop
    @assert length(candidates) > 0
    patch_p = frame[p-rₚ:p+rₚ] # pre-allocation into a contiguous memory layout
    dist = Vector{Float32}(undef, length(candidates))

    for p in R_frame
        patch_p .= @view frame[p-rₚ:p+rₚ]

        q_start = max(first(R_ref), p - S.search_radius)
        q_stop = min(last(R_ref), p + S.search_radius)
        candidates = q_start:S.search_stride:q_stop
        n = length(candidates)
        for k in 1:n
            q = candidates[k]
            patch_q = @view ref[q-rₚ:q+rₚ]
            dist[k] = S.f(patch_p, patch_q)
        end
        matches[p-rₚ, :] .= @view candidates[partialsortperm(view(dist, 1:n), 1:num_patches)]
    end

    if offset
        matches .-= R_frame
    end
    return OffsetArray(matches, (rₚ.I..., 0))
end


text:
module BondGraphs

import LightGraphs as lg
import Base: RefValue, eltype, show, in, ==, getproperty

using StaticArrays
using ModelingToolkit
using DifferentialEquations
using SymbolicUtils
using SymbolicUtils.Rewriters
using DataStructures
using Setfield
using Catalyst

export AbstractNode, Component, Junction, Port, Bond, BondGraph, BondGraphNode,
EqualEffort, EqualFlow,

type, name, vertex, set_vertex!, freeports, numports, bondgraph, nodes, bonds,
srcnode, dstnode, getnodes, getbonds,

new, add_node!, remove_node!, connect!, disconnect!,
swap!, insert_node!, merge_nodes!, simplify_junctions!,

cr, params, state_vars, set_param!, set_initial_value!,
default_value, equations, simulate

include("basetypes.jl")
#include("basetypefunctions.jl")
include("graphfunctions.jl")
include("construction.jl")
include("components.jl")
include("equations.jl")
include("conversion.jl")

end

max_stars_repo_name: jedforrest/BondGraphs

max_issues_repo_name: jedforrest/BondGraphs

max_forks_repo_name: jedforrest/BondGraphs


text:
function reproduce!(s::Simulation, μ)
	#indexarray = []
	nsteps = s.steptotal[1]
	r = .1
	for i in eachindex(s.cells)
		if s.cells[i].d == nsteps
			new_x1 = s.cells[i].p[1] + r*cosd(s.cells[i].o)
			new_y1 = s.cells[i].p[2] + r*sind(s.cells[i].o)
			new_x2 = s.cells[i].p[1] - r*cosd(s.cells[i].o)
			new_y2 = s.cells[i].p[2] - r*sind(s.cells[i].o)

			t_division1 = rand(DiscreteUniform(100+nsteps,μ+nsteps),1)
			t_division2 = rand(DiscreteUniform(100+nsteps,μ+nsteps),1)

			newcell1 = Cell([new_x1,new_y1],s.cells[i].v,t_division1[1],s.cells[i].o,birth=nsteps,motherID=i,sisterID=length(s.cells)+2)
			newcell2 = Cell([new_x2,new_y2],-s.cells[i].v,t_division2[1],s.cells[i].o,birth=nsteps,motherID=i,sisterID=length(s.cells)+1)

			# add the cell that has divided to the array of indices to be removed
			#push!(indexarray, i)
			# add the new cells
			push!(s.cells, newcell1)
			push!(s.cells, newcell2)

			s.cells[i].daughterID[:] .= [length(s.cells)+1,length(s.cells)+2]
		else
			nothing
		end
	end
	# remove the divided cells from array of cells
	#deleteat!(s.cells, indexarray)
end

max_stars_repo_name: andi-bergeson/simgym

max_issues_repo_name: andi-bergeson/simgym

max_forks_repo_name: andi-bergeson/simgym


text:
# *****************************************************************************
# Written by Ritchie Lee, ritchie.lee@sv.cmu.edu
# *****************************************************************************
# Copyright ã 2015, United States Government, as represented by the
# Administrator of the National Aeronautics and Space Administration. All
# rights reserved.  The Reinforcement Learning Encounter Simulator (RLES)
# platform is licensed under the Apache License, Version 2.0 (the "License");
# you may not use this file except in compliance with the License. You
# may obtain a copy of the License at
# http://www.apache.org/licenses/LICENSE-2.0. Unless required by applicable
# law or agreed to in writing, software distributed under the License is
# distributed on an "AS IS" BASIS, WITHOUT WARRANTIES OR CONDITIONS OF ANY
# KIND, either express or implied. See the License for the specific language
# governing permissions and limitations under the License.
# _____________________________________________________________________________
# Reinforcement Learning Encounter Simulator (RLES) includes the following
# third party software. The SISLES.jl package is licensed under the MIT Expat
# License: Copyright (c) 2014: Youngjun Kim.
# Permission is hereby granted, free of charge, to any person obtaining a copy
# of this software and associated documentation files (the "Software"), to
# deal in the Software without restriction, including without limitation the
# rights to use, copy, modify, merge, publish, distribute, sublicense, and/or
# sell copies of the Software, and to permit persons to whom the Software is
# furnished to do so, subject to the following conditions:
# The above copyright notice and this permission notice shall be included in
# all copies or substantial portions of the Software. THE SOFTWARE IS PROVIDED
# "AS IS", WITHOUT WARRANTY OF ANY KIND, EXPRESS OR IMPLIED, INCLUDING BUT
# NOT LIMITED TO THE WARRANTIES OF MERCHANTABILITY, FITNESS FOR A PARTICULAR
# PURPOSE AND NONINFRINGEMENT. IN NO EVENT SHALL THE AUTHORS OR COPYRIGHT
# HOLDERS BE LIABLE FOR ANY CLAIM, DAMAGES OR OTHER LIABILITY, WHETHER IN AN
# ACTION OF CONTRACT, TORT OR OTHERWISE, ARISING FROM, OUT OF OR IN
# CONNECTION WITH THE SOFTWARE OR THE USE OR OTHER DEALINGS IN THE SOFTWARE.
# *****************************************************************************

"""
Draw random samples and compare get_fitness() on full and pruned versions.
Every update to bestfit should be simultaneous. Produce counterexample if they differ.
Main entry: acasx_prunecheck()
"""
module ACASX_PruneCheck

export setup, runtest, acasx_prunecheck

using ACASXProblem, DerivationTrees

function setup(;data::AbstractString="libcas098small_10K",
               maxsteps::Int64=20)

  problem = ACASXClustering(data)
  grammar = get_grammar(problem)
  tree_params = DerivTreeParams(grammar, maxsteps)
  tree = DerivationTree(tree_params)

  return problem, tree
end

function runtest{T}(problem::ACASXClustering{T}, tree::DerivationTree,
                    n_samples::Int64=50000, earlystop_div::Int64=10)
  bestfit = realmax(Float64)
  for i = 1:n_samples
    rand!(tree)
    expr = get_expr(tree)

    #side-by-side comparison
    fullfit = get_fitness(problem, expr)
    prunefit = get_fitness(problem, expr, bestfit, realmax(Float64), earlystop_div)

    println("iter=$i of $(n_samples), bestfit=$(signif(bestfit, 4)), fullfit=$(signif(fullfit, 4)), prunefit=$(signif(prunefit, 4)), expr=$(string(expr))")
    if fullfit < bestfit
      bestfit = fullfit

      #on updates, these should be the same. If not, return counterexample
      if prunefit != fullfit
        return tree, expr, fullfit, prunefit
      end
    end
  end
  return tree, nothing, -1.0, -1.0
end

function acasx_prunecheck(; seed=1)
  srand(seed)
  problem, tree = setup()
  tree, expr, fullfit, prunefit = runtest(problem, tree, 50000, 10)
  @show fullfit
  @show prunefit
  @show expr
  return tree, expr, fullfit, prunefit
end

end #module


text:
"""
a simplified Boyer Moore search implementation

this is basically the boyer moore horspool algorithm, which restricts itself
to the "bad character" heuristics.

If the last byte matches, the skip is taken as maximum of the
"bad char" skip and the "good char" skip like in the bloom algorithms.
From a boyer moore viewpoint, this is a simplified special case of the
"good prefix" heuristics with a 1-byte-prefix.

There are more sophisticated implementations which care for worst case
performance and a better "good character" skip.

This implementation is not meant to earn a performance crown, it shall indicate
if (and at which pattern length) bloom filter search is comparable to
boyer moore search.
"""
function bmsearch end

function bmsearch(t::SearchSequence)
    size = sizeof(t) # no preprocessing restriction on huge patterns
    skip = size # good char skip like in bloom
    tlast = _nthbyte(t,size)
    badchar = zeros(Int,256) # for small patterns, UInt8 would be faster
    for j in 1:size-1
        badchar[_nthbyte(t,j)+1] = j
        if _nthbyte(t,j)==tlast
            skip = size-j
        end
    end
    badchar[tlast+1] = size
    return (t,badchar,skip,tlast)
end

function bmsearch(s::SearchSequence, t::SearchSequence, i::Integer)
    bmsearch(s,bmsearch(t),i)
end

function bmsearch(s::SearchSequence,p::Tuple, i::Integer,sv::MaybeVector=nothing)
    (t,badchar,skip,tlast) = p
    DOSTATS = !(sv isa Nothing)
    n = sizeof(t)
    m = sizeof(s)

    if n == 0
        return 1 <= i <= m+1 ? max(1, i) : 0
    elseif m == 0
        return 0
    elseif n == 1
        return something(findnext(isequal(_nthbyte(t,1)), s, i), 0)
    end

    w = m - n
    if w < 0 || i - 1 > w
        return 0
    end

    if DOSTATS loops = 0 end
    if DOSTATS bloomtests = 0 end # here: badchar skip > simple goodsuffix skip
    if DOSTATS bloomskips = 0 end # here: badchar skips >1
    i += n-1
    while i <= m
        if DOSTATS loops += 1 end
        if (silast = _nthbyte(s,i)) != tlast
            i += n - badchar[silast+1] # guaranteed >0
            if DOSTATS && badchar[silast+1]>1; bloomskips += 1 end
        else
            # check candidate
            j = n-1
            while (sij=_nthbyte(s,i-n+j) == _nthbyte(t,j))
                # match found?
                if j == 1
                    if DOSTATS recordcase(sv, loops, bloomtests, bloomskips, 0, skip) end
                    return i-n+1
                end
                j -= 1
            end
            # mismatch at j
            i += max(skip,j - badchar[sij+1])
            if DOSTATS && j - badchar[sij+1]>skip; bloomtests+=1;end
        end
    end
    if DOSTATS recordcase(sv, loops, bloomtests, bloomskips, 0, skip) end
    0
end


text:
#=============================#
#== put needed modules here ==#
#=============================#
import JuliaChem

#================================#
#== JuliaChem execution script ==#
#================================#
Base.@ccallable function julia_main(ARGS::Vector{String})::Cint
  #== initialize JuliaChem runtime ==#
  JuliaChem.initialize()

  #== read in input file ==#
  input_file = ARGS[1]
  molecule, driver, model, keywords = JuliaChem.JCInput.run(input_file)

  #== generate basis set ==#
  basis = JuliaChem.JCBasis.run(molecule, model)

  #== perform scf calculation ==#
  if (driver == "energy")
    if (model["method"] == "RHF")
      scf = JuliaChem.JCRHF.run(basis, molecule, keywords)
    end
  end

  #== finalize JuliaChem runtime ==#
  JuliaChem.finalize()

  #== we are done ==#
  return 0
end


text:
include("preamble.jl")

@includetests "LyceumCore.jl"


text:
# This file was generated, do not modify it. # hide
rfr = @load RandomForestRegressor pkg=DecisionTree

rfr_m = machine(rfr, X, y);


text:
function FESCHatEpsilonReader(arg0::JString, arg1::jdouble, arg2::jdouble, arg3::OceanLoadDeformationCoefficients, arg4::Map)
    return FESCHatEpsilonReader((JString, jdouble, jdouble, OceanLoadDeformationCoefficients, Map), arg0, arg1, arg2, arg3, arg4)
end

function can_add(obj::OceanTidesReader, arg0::jint, arg1::jint)
    return jcall(obj, "canAdd", jboolean, (jint, jint), arg0, arg1)
end

function get_max_parse_degree(obj::OceanTidesReader)
    return jcall(obj, "getMaxParseDegree", jint, ())
end

function get_max_parse_order(obj::OceanTidesReader)
    return jcall(obj, "getMaxParseOrder", jint, ())
end

function get_supported_names(obj::OceanTidesReader)
    return jcall(obj, "getSupportedNames", JString, ())
end

function get_waves(obj::OceanTidesReader)
    return jcall(obj, "getWaves", List, ())
end

function load_data(obj::FESCHatEpsilonReader, arg0::InputStream, arg1::JString)
    return jcall(obj, "loadData", void, (InputStream, JString), arg0, arg1)
end

function set_max_parse_degree(obj::OceanTidesReader, arg0::jint)
    return jcall(obj, "setMaxParseDegree", void, (jint,), arg0)
end

function set_max_parse_order(obj::OceanTidesReader, arg0::jint)
    return jcall(obj, "setMaxParseOrder", void, (jint,), arg0)
end

function still_accepts_data(obj::OceanTidesReader)
    return jcall(obj, "stillAcceptsData", jboolean, ())
end


text:
using PyPlot, ElasticWave

# test add born source and its adjoint
path_pv  = "/Users/wenlei/Desktop/pv.bin"
path_spt = "/Users/wenlei/Desktop/spt.bin"
function ConvertBornSource2Spts(path_spt::String, path_pv::String, dm::Array{Float64,2}, du::Array{Float64,2})
    (nz, nx, ext, iflag, dt, ns) = InfoPv(path_pv)
    if iflag == 1
       Nz = nz + 2*ext
       upper = ext
     elseif iflag == 2
       Nz = nz +   ext
       upper = 0
    end
    Nx = nx + 2*ext
    # add source to Txx
    fid = open(path_spt, "w");
    write(fid, nz, nx, ext, iflag, dt)
    for it = 1 : ns
        (vxx, vxz, vzx, vzz) = ReadPv(path_pv, it)
        Txx = zeros(Nz, Nx)
        Txx[upper+1:nz+upper, ext+1:nx+ext] = 1/2*((vxx+vzz).*dm + 2*vxx.*du)
        Txx = vec(Txx)
        Tzz = zeros(Nz, Nx)
        Tzz[upper+1:nz+upper, ext+1:nx+ext] = 1/2*((vxx+vzz).*dm + 2*vzz.*du)
        Tzz = vec(Tzz)
        Txz = zeros(Nz, Nx)
        Txz[upper+1:nz+upper, ext+1:nx+ext] = 1/2*((vxz+vzx).*du)
        Txz = vec(Txz)
        write(fid, zeros(Nz*Nx))
        write(fid, zeros(Nz*Nx))
        write(fid, zeros(Nz*Nx))
        write(fid, zeros(Nz*Nx))
        write(fid, Txx         )
        write(fid, Txx         )
        write(fid, Tzz         )
        write(fid, Tzz         )
        write(fid, Txz         )
        write(fid, Txz         )
    end
    close(fid)
    return nothing
end

dm = 100*randn(nz, nx); du = 100*randn(nz, nx);
ConvertBornSource2Spts(path_spt, path_pv, dm, du);

function makeRandSnapShots(path_spt, nz, nx, ext, iflag, dt, nt)
    if iflag == 1
       Nz = nz + 2*ext
       upper = ext
     elseif iflag == 2
       Nz = nz +   ext
       upper = 0
    end
    Nx = nx + 2*ext
    fid = open(path_spt, "w")
    write(fid, nz, nx, ext, iflag, dt);
    for it = 1 : ns
        write(fid, randn(Nz*Nx*10))
    end
    close(fid)
    return nothing
end
(nz, nx, ext, iflag, dt, ns) = InfoPv(path_pv)
path_spt1 = "/Users/wenlei/Desktop/spt1.bin"
makeRandSnapShots(path_spt1, nz, nx, ext, iflag, dt, ns)

(dm1, du1) = spts2Born(path_spt1, path_pv)
tmp = dot(vec(dm), vec(dm1)) + dot(vec(du), vec(du1))
tmp1 = IptSpts(path_spt, path_spt1)

# =======================================================

nz = 233; nx = 277; ext= 30;   iflag = 2;
dx = 5. ; dz = 5. ; dt = 5e-4; tmax = 0.5; f0=30.;
vp =4000.*ones(nz,nx);  vs=2000.*ones(nz,nx); rho=2.5*ones(nz,nx);
fidMtx = CreateFidMtx(nz, nx, ext, iflag, vp, vs, rho, dz, dx, dt, f0);

path_pv  = "/Users/wenlei/Desktop/pv.bin"
path_spt = "/Users/wenlei/Desktop/spt.bin"
dm = 100*randn(nz, nx); du = 100*randn(nz, nx);
ConvertBornSource2Spts(path_spt, path_pv, dm, du);
irx = collect(1:5:nx); irz = 2*ones(Int64, length(irx));
shotv = MultiStepForward(irz, irx, path_spt, fidMtx)
imshow(shotv.Vx, aspect=0.04)

nt = shotv.nt
shotv_adj = InitShotV(0, 0, nz, nx, ext, iflag, irz, irx, 0.0, dt, nt)
wlet = Ricker(f0, dt); nw = length(wlet)
for ir = 1 : length(irz)
    if rand() < 0.6
       tmp = randn() * wlet
       inds = floor(Int64, rand()*(nt-nw))
       shotv_adj.Vx[inds:inds+nw-1 ,ir] = tmp[:]
    end
end
imshow(shotv_adj.Vx, aspect=0.04);

path_spt_adj = "/Users/wenlei/Desktop/spt_adj.bin"
MultiStepAdjoint(path_spt_adj, shotv_adj, fidMtx)
ReverseOrderSnapShots(path_spt_adj)

tmp1 = IptSpts(path_spt_adj, path_spt)
tmp  = dot(vec(shotv.Vx), vec(shotv_adj.Vx)) + dot(vec(shotv.Vz), vec(shotv_adj.Vz))

max_stars_repo_name: WneleiGao/ElasticWave


text:
module MPIClusterManagers

export MPIManager, launch, manage, kill, procs, connect, mpiprocs, @mpi_do, TransportMode, MPI_ON_WORKERS, TCP_TRANSPORT_ALL, MPI_TRANSPORT_ALL, MPIWorkerManager

using Distributed, Serialization
import MPI

import Base: kill
import Sockets: connect, listenany, accept, IPv4, getsockname, getaddrinfo, wait_connected, IPAddr

include("workermanager.jl")
include("mpimanager.jl")
include("worker.jl")
include("mpido.jl")

end # module


text:
################################################################################
#  Copyright 2021, Marta Vanin                                                 #
################################################################################
function run_case_study_C(path_to_result_csv, nlsolver, linsolver; set_rescaler = 100, power_base::Float64=100.0)

    # Input data
    rm_transfo = true
    rd_lines = true

    season = "summer"
    time = 144
    elm = ["load", "pv"]
    pfs = [0.95, 0.90]

    ################################################################################

    # Set path
    msr_path = joinpath(mktempdir(),"temp.csv")

    # Set solve
    solver = _PMD.optimizer_with_attributes(nlsolver...)

    lin_se_solver = _PMD.optimizer_with_attributes(linsolver...)# _PMD.optimizer_with_attributes(Gurobi.Optimizer, "TimeLimit"=>180.0)

    df = _DF.DataFrame(ntw=Int64[], fdr=Int64[], solve_time=Float64[], n_bus=Int64[],
    termination_status=String[], objective=Float64[], criterion=String[], rescaler = Float64[], eq_model = String[], pbase=Float64[])

    for ntw in 1:25 for fdr in 1:10
        data_path = _PMDSE.get_enwl_dss_path(ntw, fdr)
        if !isdir(dirname(data_path)) break end

        # Load the data
        data = _PMD.parse_file(_PMDSE.get_enwl_dss_path(ntw, fdr),data_model=_PMD.ENGINEERING);
        if rm_transfo _PMDSE.rm_enwl_transformer!(data) end
        if rd_lines _PMDSE.reduce_enwl_lines_eng!(data) end
        data["settings"]["sbase_default"] = power_base

        # Insert the load profiles
        _PMDSE.insert_profiles!(data, season, elm, pfs, t = time)

        # Transform data model
        data = _PMD.transform_data_model(data);

        # Solve the ivr  power flow for linear
        pf_results_ivr = _PMD.solve_mc_pf(data, _PMD.IVRUPowerModel, solver)
        _PMDSE.write_measurements!(_PMD.IVRUPowerModel, data, pf_results_ivr, msr_path)
        _PMDSE.add_measurements!(data, msr_path, actual_meas = false, seed = 1)
        _PMDSE.assign_start_to_variables!(data)
        _PMDSE.update_all_bounds!(data; v_min = 0.8, v_max = 1.2, pg_min=-1.0, pg_max = 1.0, qg_min=-1.0, qg_max=1.0, pd_min=-1.0, pd_max=1.0, qd_min=-1.0, qd_max=1.0 )

        # Set se settings
        data["se_settings"] = Dict{String,Any}("criterion" => "rwlav",
                                "rescaler" => set_rescaler)

        linear_se_results = _SEF.run_linear_ivr_red_mc_se(data, lin_se_solver)

        pf_results_nl = _PMD.solve_mc_pf(data, _PMD.ACPUPowerModel, solver)
        _PMDSE.write_measurements!(_PMD.ACPUPowerModel, data, pf_results_nl, msr_path)
        _PMDSE.add_measurements!(data, msr_path, actual_meas = false, seed = 1)
        _PMDSE.assign_start_to_variables!(data)

        # Set se settings
        data["se_settings"] = Dict{String,Any}("criterion" => "rwlav",
                                    "rescaler" => set_rescaler)

        se_results = _PMDSE.solve_mc_se(data, _PMDSE.ReducedIVRUPowerModel, solver)

        # PRINT
        push!(df, [ntw, fdr, linear_se_results["solve_time"], length(data["bus"]),
                string(linear_se_results["termination_status"]),
                linear_se_results["objective"], "rwlav", set_rescaler, "linIVR", power_base])

        push!(df, [ntw, fdr, se_results["solve_time"], length(data["bus"]),
                string(se_results["termination_status"]),
                se_results["objective"], "rwlav", set_rescaler, "rIVR", power_base])

    end end #loop through feeder and network
    CSV.write(path_to_result_csv, df)
end