text:
abstract type AbstractDEBCore end

"""
    DEBCore(y_V_E, y_E_C, y_E_N, n_N_V, n_N_E, w_V)

Core DEB model parameters.

$(FIELDDOCTABLE)
"""
@columns struct DEBCore{MoMoD,MoMo,GMo} <: AbstractDEBCore
    # Field        | Default | Unit            | Bounds       | Log   | Description
    j_E_mai::MoMoD | 0.01    | mol*mol^-1*d^-1 | (1e-4, 1.0)  | true  | "Specific somatic maintenance costs"
    y_V_E::MoMo    | 0.7     | _               | (0.0, 1.0)   | _     | "Yield from reserve to structure"
    y_E_C::MoMo    | 0.7     | _               | (1e-6, 1.0)  | false | "Yield from C-reserve to general reserve"
    y_E_N::MoMo    | 30.0    | _               | (1.0, 50.0)  | false | "Yeild from N-reserve to general reserve"
    n_N_V::MoMo    | 0.03    | _               | (0.0, 0.1)   | _     | "Nitrogen per Carbon in structure"
    n_N_E::MoMo    | 0.025   | _               | (0.0, 0.1)   | _     | "Nitrogen per Carbon in reserve"
    w_V::GMo       | 25.0    | g*mol^-1        | (15.0, 40.0) | _     | "Mol-weight of shoot structure"
    # w_N::GMo     | 25.0    | g*mol^-1        | (15.0, 40.0) | _     | "Mol-weight of shoot N-reserve"
    # w_C::GMo     | 25.0    | g*mol^-1        | (12.0, 40.0) | _     | "Mol-weight of shoot C-reserve"
    # w_E::GMo     | 25.0    | g*mol^-1        | (15.0, 40.0) | _     | "Mol-weight of shoot reserve"
end

for fn in fieldnames(DEBCore)
    @eval $fn(p::DEBCore) = p.$fn
end

"""
    growth!(o::AbstractOrgan, u)
    growth!(p::DEBCore, o::AbstractOrgan, u)

Allocates reserves to growth flux, generalised for any number of reserves.

Where `o` is the `Organ`, and `u` is the current state parameters
"""
function growth! end
growth!(o, u) = growth!(core_pars(o), o, u)
growth!(p::DEBCore, o, u) = begin
    J = flux(o)
    J[:V,:gro] = growth = rate(o) * u[:V]
    drain = (1 / y_V_E(p)) * growth
    product = growth_production!(o, growth)
    reserve_drain!(o, :gro, drain)
end

"""
    maintenence!(o::AbstractOrgan, u)
    maintenence!(p::DEBCore, o::AbstractOrgan, u)

Allocates reserve drain due to maintenance, generalised for any number of reserves.

Maintenance is temperature dependent.

Where `o` is the `Organ`, and `u` is the current state parameters
"""
function maintenance! end
maintenence!(o, u) = maintenence!(core_pars(o), o, u)
maintenence!(p::DEBCore, o, u) = begin
    drain = j_E_mai(p) * tempcorrection(o) * u[:V]
    reserve_drain!(o, :mai, drain)
    maintenance_production!(o, u)
end

corrected_j_E_mai(o) = j_E_mai(o) * tempcorrection(o)


text:
#
println("Ce) Test of the PlasmaShift module with ASF from an internally generated multiplet.")
#
setDefaults("print summary: open", "zzz-PlasmaShift.sum")
wa = Atomic.Computation("xx",  JAC.Nuclear.Model(26.); properties=[JAC.Plasma],
                        configs=[Configuration("[Ne] 3s^2 3p^5"), Configuration("[Ne] 3s 3p^6")],
                        plasmaSettings=PlasmaShift.Settings(PlasmaShift.DebyeHueckel(), 0.25, 0., 0) )

wb = perform(wa)
setDefaults("print summary: close", "")

ext: jl


text:
using Flux3D
using Documenter

makedocs(;
    modules=[Flux3D],
    authors="Nirmal P. Suthar <nirmalps@iitk.ac.in>",
    repo="https://github.com/nirmal-suthar/Flux3D.jl/blob/{commit}{path}#L{line}",
    sitename="Flux3D.jl",
    format=Documenter.HTML(;
        prettyurls=get(ENV, "CI", "false") == "true",
        canonical="https://nirmal-suthar.github.io/Flux3D.jl",
        assets=String[],
    ),
    pages=[
        "Home" => "index.md",
    ],
)

deploydocs(;
    repo="github.com/nirmal-suthar/Flux3D.jl",
)


text:
module Testing

import ..NodalDG

using Test

@test_throws AssertionError NodalDG.Mesh1D([0, 1, 2],
                                           [[1 2]; [2 3]; [3 4]])

@test_throws AssertionError NodalDG.Mesh1D([0, 1, 2],
                                           [[0 1]; [1 2]])

let mesh = NodalDG.Mesh1D([0, 1, 2, 3],
                          [[1 2]; [2 3]; [3 4]])
    @test mesh.nvertices == 4
    @test isapprox(mesh.vertex_coordinates, [0, 1, 2, 3])
    @test mesh.nelements == 3
    @test mesh.element_vertices == [[1 2]; [2 3]; [3 4]]
end

end


text:
# *********************************************************************************
# REopt, Copyright (c) 2019-2020, Alliance for Sustainable Energy, LLC.
# All rights reserved.
#
# Redistribution and use in source and binary forms, with or without modification,
# are permitted provided that the following conditions are met:
#
# Redistributions of source code must retain the above copyright notice, this list
# of conditions and the following disclaimer.
#
# Redistributions in binary form must reproduce the above copyright notice, this
# list of conditions and the following disclaimer in the documentation and/or other
# materials provided with the distribution.
#
# Neither the name of the copyright holder nor the names of its contributors may be
# used to endorse or promote products derived from this software without specific
# prior written permission.
#
# THIS SOFTWARE IS PROVIDED BY THE COPYRIGHT HOLDERS AND CONTRIBUTORS "AS IS" AND
# ANY EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED TO, THE IMPLIED
# WARRANTIES OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR PURPOSE ARE DISCLAIMED.
# IN NO EVENT SHALL THE COPYRIGHT HOLDER OR CONTRIBUTORS BE LIABLE FOR ANY DIRECT,
# INDIRECT, INCIDENTAL, SPECIAL, EXEMPLARY, OR CONSEQUENTIAL DAMAGES (INCLUDING,
# BUT NOT LIMITED TO, PROCUREMENT OF SUBSTITUTE GOODS OR SERVICES; LOSS OF USE,
# DATA, OR PROFITS; OR BUSINESS INTERRUPTION) HOWEVER CAUSED AND ON ANY THEORY OF
# LIABILITY, WHETHER IN CONTRACT, STRICT LIABILITY, OR TORT (INCLUDING NEGLIGENCE
# OR OTHERWISE) ARISING IN ANY WAY OUT OF THE USE OF THIS SOFTWARE, EVEN IF ADVISED
# OF THE POSSIBILITY OF SUCH DAMAGE.
# *********************************************************************************
slope(x1, y1, x2, y2) = (y2 - y1) / (x2 - x1)

intercept(x1, y1, x2, y2) = y2 - slope(x1, y1, x2, y2) * x2

function insert_u_bp(xp_array_incent, yp_array_incent, region, u_xbp, u_ybp, p, u_cap)
    push!(xp_array_incent[region], u_xbp)
    push!(yp_array_incent[region], u_ybp - u_ybp * p + u_cap)
    return xp_array_incent, yp_array_incent
end

function insert_p_bp(xp_array_incent, yp_array_incent, region, p_xbp, p_ybp, u, p_cap)
    push!(xp_array_incent[region], p_xbp)
    push!(yp_array_incent[region], p_ybp - (p_cap + p_xbp * u))
    return xp_array_incent, yp_array_incent
end

function insert_u_after_p_bp(xp_array_incent, yp_array_incent, region, u_xbp, u_ybp, p, p_cap, u_cap)
    xp_array_incent[region].append(u_xbp)
    if p_cap == 0
        push!(yp_array_incent[region], u_ybp - (p * u_ybp + u_cap))
    else
        push!(yp_array_incent[region], u_ybp - (p_cap + u_cap))
    end
    return xp_array_incent, yp_array_incent
end

function insert_p_after_u_bp(xp_array_incent, yp_array_incent, region, p_xbp, p_ybp, u, u_cap, p_cap)
    push!(xp_array_incent[region], p_xbp)
    if u_cap == 0
        push!(yp_array_incent[region], p_ybp - (p_cap + u * p_xbp))
    else
        push!(yp_array_incent[region], p_ybp - (p_cap + u_cap))
    end
    return xp_array_incent, yp_array_incent
end

"""
    cost_curve(tech::AbstractTech, financial::Financial)

Determine the cost curve segments (x and y points) accounting for:
    1. tax incentives
    2. rebates
    3. (size, cost) samples

Assumes that `tech` has the following attributes:
- installed_cost_per_kw
- federal_itc_pct
- federal_rebate_per_kw
- state_ibi_pct
- state_ibi_max
- state_rebate_per_kw
- state_rebate_max
- utility_ibi_pct
- utility_ibi_max
- utility_rebate_per_kw
- utility_rebate_max
- macrs_option_years
- macrs_bonus_pct
- macrs_itc_reduction

Optional attributes of `tech` are:
- existing_kw
- tech_sizes_for_cost_curve
"""
function cost_curve(tech::AbstractTech, financial::Financial)
    big_number = 1.0e10
    T = typeof(tech)

    if nameof(T) in [:Boiler, :Elecchl]  # zero cost technologies
        cap_cost_slope = [0.0]
        cap_cost_yint = [0.0]
        cost_curve_bp_x = [0.0, big_number]
        n_segments = 1
        return cap_cost_slope, cost_curve_bp_x, cap_cost_yint, n_segments
    end

    regions = ["utility", "state", "federal", "combined"]
    cap_cost_slope = Real[]
    cost_curve_bp_x = [0.0]
    cost_curve_bp_y = [0.0]
    cap_cost_yint = Real[]
    n_segments = 1

    existing_kw = 0.0
    if :existing_kw in fieldnames(T)
        existing_kw = tech.existing_kw
    end

    tech_incentives = Dict(
        "federal" => Dict(),
        "state" => Dict(),
        "utility" => Dict(),
    )
    if nameof(T) in [:Generator, :Absorpchl]  # for generator and absorption chiller, there are no incentives
        for region in keys(tech_incentives)
            tech_incentives[region]["%"] = 0.0
            tech_incentives[region]["%_max"] = 0.0
            tech_incentives[region]["rebate"] = 0.0
            tech_incentives[region]["rebate_max"] = 0.0
        end
    else
        # NOTE REopt incentive calculation works best if "unlimited" incentives are entered as 0
        tech_incentives["federal"]["%"] = tech.federal_itc_pct
        tech_incentives["federal"]["%_max"] = 0
        tech_incentives["state"]["%"] = tech.state_ibi_pct
        tech_incentives["state"]["%_max"] = tech.state_ibi_max ==  big_number ? 0 : tech.state_ibi_max
        tech_incentives["utility"]["%"] = tech.utility_ibi_pct
        tech_incentives["utility"]["%_max"] = tech.utility_ibi_max == big_number ? 0 : tech.utility_ibi_max

        tech_incentives["federal"]["rebate"] = tech.federal_rebate_per_kw
        tech_incentives["federal"]["rebate_max"] = 0
        tech_incentives["state"]["rebate"] = tech.state_rebate_per_kw
        tech_incentives["state"]["rebate_max"] = tech.state_rebate_max ==  big_number ? 0 : tech.state_rebate_max
        tech_incentives["utility"]["rebate"] = tech.utility_rebate_per_kw
        tech_incentives["utility"]["rebate_max"] = tech.utility_rebate_max  ==  big_number ? 0 : tech.utility_rebate_max
    end

    # Intermediate Cost curve
    xp_array_incent = Dict("utility" => [0.0, big_number])
    yp_array_incent = Dict("utility" => [0.0, big_number * tech.installed_cost_per_kw])  # [$]

    # New input of tech_sizes_for_cost_curve to be associated with tech.installed_cost_per_kw with same type and length
    if :tech_sizes_for_cost_curve in fieldnames(T)
        if !isempty(tech.tech_sizes_for_cost_curve)
            if length(tech.tech_sizes_for_cost_curve) == 1
                yp_array_incent["utility"] = [0.0, big_number * tech.installed_cost_per_kw[1]]  # [$]
            else  # length(tech.tech_sizes_for_cost_curve) > 1
                xp_array_incent["utility"] = []
                if tech.tech_sizes_for_cost_curve[1] != 0 # Append a 0 to the front of the list if not included (we"ll assume that it has a 0 y-intercept below)
                    push!(xp_array_incent["utility"], 0)
                end
                push!(xp_array_incent["utility"], tech.tech_sizes_for_cost_curve...)  # [$]  # Append list of sizes for cost curve [kW]
                push!(xp_array_incent["utility"], big_number)  # Append big number size to assume same cost as last input point

                if tech.tech_sizes_for_cost_curve[1] == 0
                    yp_array_incent["utility"] = [tech.installed_cost_per_kw[1], tech.tech_sizes_for_cost_curve[2:end] .* tech.installed_cost_per_kw[2:end]...]  # [$]
                    # tech.installed_cost_per_kw[1] is assumed to be in units of $, if there is tech.tech_sizes_for_cost_curve[1]=0 point
                else
                    yp_array_incent["utility"] = [0, tech.tech_sizes_for_cost_curve .* tech.installed_cost_per_kw...]
                end
                append!(yp_array_incent["utility"], big_number * tech.installed_cost_per_kw[end])  # Last cost assumed for big_number size
                cost_curve_bp_y = [yp_array_incent["utility"][1]]
            end
        end
    end

    for r in range(1, stop=length(regions)-1)  # apply incentives

        region = regions[r]  # regions = ["utility", "state", "federal", "combined"]
        next_region = regions[r + 1]

        # Apply incentives, initialize first value
        xp_array_incent[next_region] = [0.0]
        yp_array_incent[next_region] = [0.0]

        # percentage based incentives
        p = tech_incentives[region]["%"]
        p_cap = tech_incentives[region]["%_max"]

        # rebates, for some reason called "u" in REopt
        u = tech_incentives[region]["rebate"]
        u_cap = tech_incentives[region]["rebate_max"]

        # reset switches and break point counter
        switch_percentage = (p == 0 || p_cap == 0) ? true : false
        switch_rebate = (u == 0 || u_cap == 0) ? true : false

        # start at second point, first is always zero
        for point in range(2, stop=length(xp_array_incent[region]))

            # previous points
            xp_prev = xp_array_incent[region][point - 1]
            yp_prev = yp_array_incent[region][point - 1]

            # current, unadjusted points
            xp = xp_array_incent[region][point]
            yp = yp_array_incent[region][point]

            # initialize the adjusted points on cost curve
            xa = xp
            ya = yp

            # initialize break points
            u_xbp = 0.0
            u_ybp = 0.0
            p_xbp = 0.0
            p_ybp = 0.0

            if  !switch_rebate
                u_xbp = u_cap / u
                u_ybp = slope(xp_prev, yp_prev, xp, yp) * u_xbp + intercept(xp_prev, yp_prev, xp, yp)
            end

            if !switch_percentage
                p_xbp = (p_cap / p - intercept(xp_prev, yp_prev, xp, yp)) / slope(xp_prev, yp_prev, xp, yp)
                p_ybp = p_cap / p
            end

            if ((p * yp) < p_cap || p_cap == 0) && ((u * xp) < u_cap || u_cap == 0)
                ya = yp - (p * yp + u * xp)
            elseif (p * yp) < p_cap && (u * xp) >= u_cap
                if !switch_rebate
                    if u * xp != u_cap
                        xp_array_incent, yp_array_incent =
                            insert_u_bp(xp_array_incent, yp_array_incent, next_region, u_xbp, u_ybp, p, u_cap)
                    end
                    switch_rebate = true
                end
                ya = yp - (p * yp + u_cap)
            elseif (p * yp) >= p_cap && (u * xp) < u_cap
                if !switch_percentage
                    if p * yp != p_cap
                        xp_array_incent, yp_array_incent =
                            insert_p_bp(xp_array_incent, yp_array_incent, next_region, p_xbp, p_ybp, u, p_cap)
                    end
                    switch_percentage = true
                end
                ya = yp - (p_cap + xp * u)
            elseif p * yp >= p_cap && u * xp >= u_cap
                if !switch_rebate && !switch_percentage
                    if p_xbp == u_xbp
                        xp_array_incent, yp_array_incent =
                            insert_u_bp(xp_array_incent, yp_array_incent, next_region, u_xbp, u_ybp, p, u_cap)
                        switch_percentage = true
                        switch_rebate = true
                    elseif p_xbp < u_xbp
                        if p * yp != p_cap
                            xp_array_incent, yp_array_incent =
                                insert_p_bp(xp_array_incent, yp_array_incent, next_region, p_xbp, p_ybp, u,
                                            p_cap)
                        end
                        switch_percentage = true
                        if u * xp != u_cap
                            xp_array_incent, yp_array_incent =
                                insert_u_after_p_bp(xp_array_incent, yp_array_incent, next_region, u_xbp, u_ybp,
                                                    p, p_cap, u_cap)
                        end
                        switch_rebate = true
                    else
                        if u * xp != u_cap
                            xp_array_incent, yp_array_incent =
                                insert_u_bp(xp_array_incent, yp_array_incent, next_region, u_xbp, u_ybp, p, u_cap)
                        end
                        switch_rebate = true
                        if p * yp != p_cap
                            xp_array_incent, yp_array_incent =
                                insert_p_after_u_bp(xp_array_incent, yp_array_incent, next_region, p_xbp, p_ybp,
                                                    u, u_cap, p_cap)
                        end
                        switch_percentage = true
                    end
                elseif switch_rebate && !switch_percentage
                    if p * yp != p_cap
                        xp_array_incent, yp_array_incent =
                            insert_p_after_u_bp(xp_array_incent, yp_array_incent, next_region, p_xbp, p_ybp, u,
                                                u_cap, p_cap)
                    end
                    switch_percentage = true
                elseif !switch_rebate && switch_percentage
                    if u * xp != u_cap
                        xp_array_incent, yp_array_incent =
                            insert_u_after_p_bp(xp_array_incent, yp_array_incent, next_region, u_xbp, u_ybp, p,
                                                p_cap, u_cap)
                    end
                    switch_rebate = true
                end

                # Finally compute adjusted values
                if p_cap == 0
                    ya = yp - (p * yp + u_cap)
                elseif u_cap == 0
                    ya = yp - (p_cap + u * xp)
                else
                    ya = yp - (p_cap + u_cap)
                end
            end

            push!(xp_array_incent[next_region], xa)
            push!(yp_array_incent[next_region], ya)

            # compute cost curve, funky logic in REopt ignores everything except xa, ya
            if region == "federal"
                push!(cost_curve_bp_x, xa)
                push!(cost_curve_bp_y, ya)
            end
        end
    end

    for seg in range(2, stop=length(cost_curve_bp_x))
        tmp_slope = round((cost_curve_bp_y[seg] - cost_curve_bp_y[seg - 1]) /
                        (cost_curve_bp_x[seg] - cost_curve_bp_x[seg - 1]), digits=0)
        tmp_y_int = round(cost_curve_bp_y[seg] - tmp_slope * cost_curve_bp_x[seg], digits=0)

        push!(cap_cost_slope, tmp_slope)
        push!(cap_cost_yint, tmp_y_int)
    end

    n_segments = length(cap_cost_slope)

    # Following logic modifies the cap cost segments to account for the tax benefits of the ITC and MACRs
    updated_cap_cost_slope = Real[]
    updated_y_intercept = Real[]

    for s in range(1, stop=n_segments)
        if cost_curve_bp_x[s + 1] <= 0
            # Not sure how else to handle this case, perhaps there is a better way to handle it?
            @error "Invalid cost curve for {$nameof(T)}. Value at index {$s} ({$cost_curve_bp_x[s + 1]}) cannot be less than or equal to 0."
        end

        # Remove federal incentives for ITC basis and tax benefit calculations
        itc = tech.federal_itc_pct
        rebate_federal = tech.federal_rebate_per_kw
        if itc == 1
            itc_unit_basis = 0
        else
            itc_unit_basis = (cap_cost_slope[s] + rebate_federal) / (1 - itc)
        end

        macrs_schedule = [0.0]
        macrs_bonus_pct = 0.0
        macrs_itc_reduction = 0.0

        if tech.macrs_option_years != 0
            macrs_bonus_pct = tech.macrs_bonus_pct
            macrs_itc_reduction = tech.macrs_itc_reduction
        end
        if tech.macrs_option_years == 5
            macrs_schedule = financial.macrs_five_year
        end
        if tech.macrs_option_years == 7
            macrs_schedule = financial.macrs_seven_year
        end
        updated_slope = effective_cost(;
            itc_basis=itc_unit_basis,  # input tech cost with incentives, but no ITC
            replacement_cost=0.0,
            replacement_year=financial.analysis_years,
            discount_rate=financial.owner_discount_pct,
            tax_rate=financial.owner_tax_pct,
            itc=itc,
            macrs_schedule = macrs_schedule,
            macrs_bonus_pct=macrs_bonus_pct,
            macrs_itc_reduction = macrs_itc_reduction,
            rebate_per_kw = rebate_federal
        )
        # The way REopt incentives currently work, the federal rebate is the only incentive that doesn't reduce ITC basis
        push!(updated_cap_cost_slope, updated_slope)
    end

    for p in range(2, stop=n_segments + 1)
        cost_curve_bp_y[p] = cost_curve_bp_y[p - 1] + updated_cap_cost_slope[p - 1] *
                                (cost_curve_bp_x[p] - cost_curve_bp_x[p - 1])
        push!(updated_y_intercept, cost_curve_bp_y[p] - updated_cap_cost_slope[p - 1] * cost_curve_bp_x[p])
    end
    cap_cost_slope = updated_cap_cost_slope
    cap_cost_yint = updated_y_intercept

    return cap_cost_slope, cost_curve_bp_x, cap_cost_yint, n_segments
end


text:
@testset "End of Week Crossover" begin
    t = 0
    @test t == (PositionVelocityTime.check_crossover(deepcopy(t)))

    t = 350000
    @test t - 604800 == PositionVelocityTime.check_crossover(deepcopy(t))

    t = -350000
    @test t + 604800 == PositionVelocityTime.check_crossover(deepcopy(t))

end

correct_raw_times =  [
    208863.7577188781,
    208863.76169284966,
    208863.76290867585,
    208863.75739086428,
    208863.7704948417]

@testset "Calc uncorrected Time" begin
    for i in 1:length(satellite_states)
        out = PositionVelocityTime.calc_uncorrected_time(satellite_states[i])
        @test out == correct_raw_times[i]
    end
end


text:
@testset "parse fenced" begin
    s = """
        A
        ```julia
        using Random
        Random.seed!(55) # hide
        a = randn()
        ```
        B
        """
    css = SubString(s, 3, 63)
    lang, rpath, code = F.parse_fenced_block(css)
    @test lang == "julia"
    @test rpath === nothing
    @test code ==
        """using Random
        Random.seed!(55) # hide
        a = randn()"""
    s = """
        A
        ```julia:ex1
        using Random
        Random.seed!(55) # hide
        a = randn()
        ```
        B
        """
    css = SubString(s, 3, 67)
    lang, rpath, code = F.parse_fenced_block(css)
    @test lang == "julia"
    @test rpath == "ex1"
    @test code ==
        """using Random
        Random.seed!(55) # hide
        a = randn()"""

    # more cases to test REGEX_CODE
    rx3 = F.CODE_3_PAT
    rx5 = F.CODE_5_PAT

    s = "```julia:ex1 A```"
    m = match(rx3, s)
    @test m.captures[2] == ":ex1"
    s = "```julia:ex1_b2 A```"
    m = match(rx3, s)
    @test m.captures[2] == ":ex1_b2"
    s = "```julia:ex1_b2-33 A```"
    m = match(rx3, s)
    @test m.captures[2] == ":ex1_b2-33"
    s = "```julia:./ex1/v1_b2 A```"
    m = match(rx3, s)
    @test m.captures[2] == ":./ex1/v1_b2"
    s = "```julia:./ex1/v1_b2.f99 A```"
    m = match(rx3, s)
    @test m.captures[2] == ":./ex1/v1_b2.f99"

    s = "`````julia:./ex1/v1_b2.f99 A`````"
    m = match(rx5, s)
    @test m.captures[2] == ":./ex1/v1_b2.f99"
end

@testset "resolve code" begin
    # no eval
    c = SubString(
        """```julia
        a = 5
        b = 7
        ```""")
    @test F.resolve_code_block(c) ==
        """
$(F.htmlesc("""a = 5
        b = 7"""))
"""
    # not julia code
    c = SubString(
        """```python:ex
        a = 5
        b = 7
        ```""")

    s = @capture_out F.resolve_code_block(c)
    @test occursin("Evaluation of non-Julia code", s)

    # should eval
    bak = pwd()
    tmp = mktempdir()
    begin
        cd(tmp)
        set_curpath("index.md")
        F.FOLDER_PATH[] = tmp
        F.def_LOCAL_VARS!()
        F.set_paths!()
        c = SubString(
            """```julia:ex
            a = 5
            b = 7
            ```""")
        r = F.resolve_code_block(c)
        @test r ==
            """
$(F.htmlesc("""a = 5
            b = 7"""))
"""
    end
    cd(bak)
end


text:
# ------------------------------------------------------------------- #
# Epstein Zin
# ------------------------------------------------------------------- #
immutable EpsteinZinAgent <: AbstractConsumer
    ρ::Float64
    α::Float64
    β::Float64
end

_unpack(ez::EpsteinZinAgent) = (ez.ρ, ez.α, ez.β)

is_recursive(a::EpsteinZinAgent) = a.ρ == a.α ? false : true

utility(a::EpsteinZinAgent, flow, continuation) =
    ((1 - a.β).*flow.^a.ρ + a.β * continuation.^a.ρ).^(1./a.ρ)

# ------------------------------------------------------------------- #
# Power Utility
# ------------------------------------------------------------------- #
immutable CRRAAgent <: AbstractConsumer
    γ::Float64
    β::Float64
end

typealias PowerUtilityAgent CRRAAgent

is_recursive(::CRRAAgent) = false
utility(a::CRRAAgent, c) = c.^(1-a.γ)./(1 - a.γ)

max_stars_repo_name: NYUEcon/GrowthModels

max_issues_repo_name: NYUEcon/GrowthModels

max_forks_repo_name: NYUEcon/GrowthModels


text:
include("network_library.jl")
#  include("roadnet_pursuer_generator_MDP.jl")

ext_rwd = 2000.
cgt_rwd = -2000.

go = original_roadnet(exit_rwd=ext_rwd,caught_rwd=cgt_rwd,sensor_rwd=-200.)
g = medium_roadnet(exit_rwd=ext_rwd,caught_rwd=cgt_rwd,sensor_rwd=-200.)
#  mdp = roadnet_with_pursuer(g,tp=0.7,d=0.9)

display_network(go,evader_locs=[1],pursuer_locs=[4],fname="original_net",scale=1.0)
display_network(g,evader_locs=[1],pursuer_locs=[4],fname="medium_net",scale=2.0)

max_stars_repo_name: COHRINT/FaMSeC

max_issues_repo_name: COHRINT/FaMSeC

max_forks_repo_name: COHRINT/FaMSeC


text:
function test9()
  sim = LTspiceSimulation("test9.asc",executablepath="")
  @test LTspice.does_circuitfilearray_file_match(sim)
  show(IOBuffer(),sim)
  compare_arrays(measurementvalues(sim),[-0.5 -1.0 1.0; 0.0 NaN NaN; 0.5 1.0 -1.0])
#  @test measurementvalues(sim)≈[-0.5 -1.0 1.0; 0.0 NaN NaN; 0.5 1.0 -1.0]
  show(IOBuffer(),sim)
end
test9()


text:
## Include needed packages
using IncrementalInference
using RoMEPlotting

# import getSample to be extended for user factor MultiModalConditional
import IncrementalInference: getSample

## create a new facor type MultiModalConditional
## FIXME, this approach is unnecessary, see `::Mixture` instead.  See Caesar.jl documentation for details.
mutable struct MultiModalConditional <: AbstractRelativeRoots
  x::Vector{Distribution}
  hypo::Categorical
  MultiModalConditional(x::Vector{<:Distribution}, p::Categorical) = new(x, p)
end
function getSample(cf::CalcFactor{<:MultiModalConditional}, N::Int=1)
  d = length(cf.factor.hypo.p)
  p = rand(cf.factor.hypo, N)
  ret = zeros(1,N)
  for i in 1:N
    ret[i] = rand(cf.factor.x[p[i]])
  end
  return (ret, p)
end

function (cf::CalcFactor{<:MultiModalConditional})(meas, x1, x2)
  #
  return meas[1] - (x2[1]-x1[1])
end

## build factor graph and populate
fg = initfg()

N=100

doors = [-20.0 0.0 20.0]
pd = kde!(doors,[2.0])
pd = resample(pd,N);
bws = getBW(pd)[:,1]
doors2 = getPoints(pd);
v1 = addVariable!(fg,:x1,ContinuousScalar,N=N)
f1 = addFactor!(fg,[v1],Prior(pd)) #, samplefnc=getSample

# not initialized
v2 = addVariable!(fg,:x2, ContinuousScalar, N=N)

mmc = MultiModalConditional([Normal(-5,0.5),Normal(5,0.5)],Categorical([0.5,0.5]))
f2 = addFactor!(fg, [:x1; :x2], mmc )

pts = approxConv(fg, :x1x2f1, :x2)

## do some plotting
meas = sampleFactor(f2,2000)
q2 = kde!(meas)
h1 = plotKDE([getBelief(v1), q2],c=["red";"green"],fill=true, xlbl="")
h2 = plotKDE(kde!(pts),fill=true,xlbl="", title="N = 100")

draw(PDF("approxconv.pdf",14cm,10cm),vstack(h1,h2))
# @async run(`evince approxconv.pdf`)


text:
module QJuliaMR

using QJuliaInterface
using QJuliaEnums
using QJuliaBlas
using QJuliaReduce
using QJuliaSolvers

using LinearAlgebra
using MPI

norm2             = QJuliaReduce.gnorm2
reDotProductNormX = QJuliaReduce.reDotProductNormX
axpyXmaz          = QJuliaBlas.axpyXmaz

verbose           = false

function solver(x::AbstractArray, b::AbstractArray, Mat::Any, MatSloppy::Any, param::QJuliaSolvers.QJuliaSolverParam_qj)

    if verbose == true; println("Running MR solver in maximum " , param.Nsteps, " step(s)."); end

    QJuliaReduce.set_blas_global_reduction(param.global_reduction)

    if MPI.Initialized() == false; error("MPI was not inititalized, copy source field to the solution."); end

    if (param.maxiter == 0) || (param.Nsteps == 0)
      if param.use_init_guess == false
        x .=@. 0.0
      end
      return
    end #if

    mixed = (param.dtype_sloppy != param.dtype)

    local r   = Vector{param.dtype}(undef, length(b))
    # now allocate sloppy fields
    local rSloppy = mixed == true ? Vector{param.dtype_sloppy}(undef, length(b)) : r
    local Ar      = typeof(rSloppy)(undef, length(rSloppy))
    #  iterated sloppy solution vector
    local xSloppy = typeof(rSloppy)(undef, length(rSloppy))

    b2 = norm2(b)  #Save norm of b
    r2 = 0.0     #if zero source then we will exit immediately doing no work

    if param.use_init_guess == true
      #r = b - Ax0 <- real
      Mat(r, x)
      r .=@. b - r
      r2 = norm2(r)
    else
      r2 = b2
      r .=@. b
      x .=@. 0.0
    end

    rSloppy .=@. r

    # if invalid residual then convergence is set by iteration count only
    stop = b2*param.tol*param.tol

    if verbose == true;println("MR: Initial residual = ", sqrt(r2));end

    converged = false; step = 0

    while converged == false

      scale = 1.0
      #set to zero sloppy solution
      xSloppy .=@. 0.0

      c2    = r2  # c2 holds the initial r2
      scale = c2 > 0.0 ? sqrt(c2) : 1.0

      # domain-wise normalization of the initial residual to prevent underflow
      if (c2 > 0.0)
	rSloppy ./=@. scale
	r2 = 1.0
      end

      k = 0
      if verbose == true;println("MR: ", step, " cycle, ",  k, " iterations, r2 = ", r2);end

      while (k < param.maxiter && r2 > 0.0)
        MatSloppy(Ar, rSloppy)

        alpha = reDotProductNormX(Ar, rSloppy)
	# x += omega*alpha*r, r -= omega*alpha*Ar, r2 = blas::norm2(r)//?
        coeff = (param.omega*alpha[1]) / alpha[2]
	axpyXmaz(coeff, rSloppy, xSloppy, Ar)
	if(verbose == true); println("MR: ", step ," cycle, ", (k+1)," iterations, <r|A|r> = ", alpha[1], "  residual norm = ", norm2(rSloppy)); end

        k += 1
      end #while k < param.maxiter && r2 > 0.0

      # Scale and sum to accumulator
      x .=@. x + scale*xSloppy

      step += 1

      if (param.compute_true_res == true || param.Nsteps > 1)
        Mat(r, x)

        r .=@. b - r
        r2 = norm2(r)

        param.true_res = sqrt(r2 / b2)

        converged = (step < param.Nsteps && r2 > stop) ? false : true

        if (param.preserve_source == false && converged == true)
          b .=@. r
        else
          rSloppy .=@. r
        end

        println("MR: ", step ," cycle, Converged after ", k , "  iterations, relative residual: true = ", sqrt(r2))

      else
        rSloppy .*= scale

        r2 = norm2(rSloppy)

        converged = (step < param.Nsteps) ? false : true

        if (param.preserve_source == false && converged == true)
          b .=@. rSloppy
        else
          r .=@. rSloppy
        end

        println("MR: ", step ," cycle, Converged after ", k , "  iterations, relative residual: true = ", sqrt(r2))

      end #if (param.compute_true_res == true || param.Nsteps > 1)

    end #while converged == false

    QJuliaReduce.reset_blas_global_reduction()

end #MR

end #QJuliaMR

max_stars_repo_name: alexstrel/QJulia

max_issues_repo_name: alexstrel/QJulia

max_forks_repo_name: alexstrel/QJulia


text:
function solve()
    n = parse(Int, readline())
    cs = readline()

    numr = count("R", cs)
    count("W", cs[1:numr])
end

println(solve())

max_stars_repo_name: aishikawa/atcoder-julia

max_issues_repo_name: aishikawa/atcoder-julia

max_forks_repo_name: aishikawa/atcoder-julia


text:
using SafeTestsets, Test, Pkg

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

function dev_subpkg(subpkg)
    subpkg_path = joinpath(dirname(@__DIR__), "lib", subpkg)
    Pkg.develop(PackageSpec(path=subpkg_path))
end

function activate_subpkg_env(subpkg)
    subpkg_path = joinpath(dirname(@__DIR__), "lib", subpkg)
    Pkg.activate(subpkg_path)
    Pkg.develop(PackageSpec(path=subpkg_path))
    Pkg.instantiate()
end

groups = if GROUP == "All"
    ["Lux", "Boltz"]
else
    [GROUP]
end

@time begin for group in groups
    @info "Testing GROUP $group"
    if group == "Lux"
        @testset "Lux.jl" begin
            @testset "Layers" begin
                @time @safetestset "Basic" begin include("layers/basic.jl") end
                @time @safetestset "Convolution" begin include("layers/conv.jl") end
                @time @safetestset "Normalization" begin include("layers/normalize.jl") end
                @time @safetestset "Recurrent" begin include("layers/recurrent.jl") end
                @time @safetestset "Dropout" begin include("layers/dropout.jl") end
            end

            @time @safetestset "Functional Operations" begin include("functional.jl") end
        end
    else
        dev_subpkg(group)
        subpkg_path = joinpath(dirname(@__DIR__), "lib", group)
        Pkg.test(PackageSpec(name=group, path=subpkg_path))
    end
end end


text:
using Jtac

i = 0
ds = record_self(MCTSPlayer(power=10), 10, game = MetaTac(), branching = 0.1, callback = () -> (global i+=1; println(i)))

max_stars_repo_name: roSievers/jtac


text:
## 14.jl
# 先頭からN行を出力

f = open("./dat/hightemp.txt")
n = parse(Int, ARGS[1])
for line in readlines(f)[1:n]
  println(chomp(line))
end


text:
cd("$(homedir())/Dropbox/Doctorate/GIT/BSyncGPGPU/")

tempFeatures <- readtable("$(homedir())/Dropbox/Doctorate/GIT/BSyncGPGPU/datasets/AllAppGPUs.csv")

df =

max_stars_repo_name: marcosamaris/BSPGPUSync

max_issues_repo_name: marcosamaris/BSPGPUSync

max_forks_repo_name: marcosamaris/BSPGPUSync


text:
export FiniteDifferenceMethod, fdm, backward_fdm, forward_fdm, central_fdm

"""
    FiniteDifferences.DEFAULT_CONDITION

The default [condition number](https://en.wikipedia.org/wiki/Condition_number) used when
computing bounds. It provides amplification of the ∞-norm when passed to the function's
derivatives.
"""
const DEFAULT_CONDITION = 100

"""
    FiniteDifferences.TINY

A tiny number added to some quantities to ensure that division by 0 does not occur.
"""
const TINY = 1e-20

forward_grid(p::Int) = 0:(p - 1)
backward_grid(p::Int) = (1 - p):0
function central_grid(p::Int)
    if isodd(p)
        return div(1 - p, 2):div(p - 1, 2)
    else
        return vcat(div(-p, 2):-1, 1:div(p, 2))
    end
end

"""
    History

A mutable type that tracks several values during adaptive bound computation.
"""
mutable struct History
    adapt::Int
    eps::Real
    bound::Real
    step::Real
    accuracy::Real

    function History(; kwargs...)
        h = new()
        for (k, v) in kwargs
            setfield!(h, k, v)
        end
        return h
    end
end

"""
    FiniteDifferenceMethod

Abstract type for all finite differencing method types.
Subtypes of `FiniteDifferenceMethod` are callable with the signature

```
method(f, x; kwargs...)
```

where the keyword arguments can be any of

* `adapt`: The number of adaptive steps to use improve the estimate of `bound`.
* `bound`: Bound on the value of the function and its derivatives at `x`.
* `condition`: The condition number. See [`DEFAULT_CONDITION`](@ref).
* `eps`: The assumed roundoff error. Defaults to `eps()` plus [`TINY`](@ref).
"""
abstract type FiniteDifferenceMethod end

function Base.show(io::IO, x::FiniteDifferenceMethod)
    @printf io "FiniteDifferenceMethod:\n"
    @printf io "  order of method:       %d\n" x.p
    @printf io "  order of derivative:   %d\n" x.q
    @printf io "  grid:                  %s\n" x.grid
    @printf io "  coefficients:          %s\n" x.coefs
    h = x.history
    if all(p->isdefined(h, p), propertynames(h))
        @printf io "  roundoff error:        %.2e\n" h.eps
        @printf io "  bounds on derivatives: %.2e\n" h.bound
        @printf io "  step size:             %.2e\n" h.step
        @printf io "  accuracy:              %.2e\n" h.accuracy
    end
end

for D in (:Forward, :Backward, :Central, :Nonstandard)
    @eval begin
        struct $D{G<:AbstractVector, C<:AbstractVector} <: FiniteDifferenceMethod
            p::Int
            q::Int
            grid::G
            coefs::C
            history::History
        end
        (d::$D)(f, x=0.0; kwargs...) = fdm(d, f, x; kwargs...)
    end
end

# The below does not apply to Nonstandard, as it has its own constructor
for D in (:Forward, :Backward, :Central)
    lcname = lowercase(String(D))
    gridf = Symbol(lcname, "_grid")
    fdmf = Symbol(lcname, "_fdm")

    @eval begin
        # Compatibility layer over the "old" API
        function $fdmf(p::Integer, q::Integer; adapt=1, kwargs...)
            _dep_kwarg(kwargs)
            return $D(p, q; adapt=adapt, kwargs...)
        end

        function $D(p::Integer, q::Integer; adapt=1, kwargs...)
            _check_p_q(p, q)
            grid = $gridf(p)
            coefs = _coefs(grid, p, q)
            hist = History(; adapt=adapt, kwargs...)
            return $D{typeof(grid), typeof(coefs)}(Int(p), Int(q), grid, coefs, hist)
        end

        @doc """
            FiniteDifferences.$($(Meta.quot(D)))(p, q; kwargs...)
            $($(Meta.quot(fdmf)))(p, q; kwargs...)

        Construct a $($lcname) finite difference method of order `p` to compute the `q`th
        derivative.
        See [`FiniteDifferenceMethod`](@ref) for more details.
        """
        ($D, $fdmf)
    end
end

"""
    FiniteDifferences.Nonstandard(grid, q; kwargs...)

An finite differencing method which is constructed based on a user-defined grid. It is
nonstandard in the sense that it represents neither forward, backward, nor central
differencing.
See [`FiniteDifferenceMethod`](@ref) for further details.
"""
function Nonstandard(grid::AbstractVector{<:Real}, q::Integer; adapt=0, kwargs...)
    p = length(grid)
    _check_p_q(p, q)
    coefs = _coefs(grid, p, q)
    hist = History(; adapt=adapt, kwargs...)
    return Nonstandard{typeof(grid), typeof(coefs)}(Int(p), Int(q), grid, coefs, hist)
end

# Check the method and derivative orders for consistency
function _check_p_q(p::Integer, q::Integer)
    q >= 0 || throw(ArgumentError("order of derivative must be nonnegative"))
    q < p || throw(ArgumentError("order of the method must be strictly greater than that " *
                                 "of the derivative"))
    # Check whether the method can be computed. We require the factorial of the
    # method order to be computable with regular `Int`s, but `factorial` will overflow
    # after 20, so 20 is the largest we can allow.
    p > 20 && throw(ArgumentError("order of the method is too large to be computed"))
    return
end

# Compute coefficients for the method
function _coefs(grid::AbstractVector{<:Real}, p::Integer, q::Integer)
    C = [g^i for i in 0:(p - 1), g in grid]
    x = zeros(Int, p)
    x[q + 1] = factorial(q)
    return C \ x
end

# Estimate the bound on the function value and its derivatives at a point
_estimate_bound(x, cond) = cond * maximum(abs, x) + TINY

"""
    fdm(m::FiniteDifferenceMethod, f, x[, Val(false)]; kwargs...) -> Real
    fdm(m::FiniteDifferenceMethod, f, x, Val(true); kwargs...) -> Tuple{FiniteDifferenceMethod, Real}

Compute the derivative of `f` at `x` using the finite differencing method `m`.
The optional `Val` argument dictates whether the method should be returned alongside the
derivative value, which can be useful for examining the step size used and other such
parameters.

The recognized keywords are:

* `adapt`: The number of adaptive steps to use improve the estimate of `bound`.
* `bound`: Bound on the value of the function and its derivatives at `x`.
* `condition`: The condition number. See [`DEFAULT_CONDITION`](@ref).
* `eps`: The assumed roundoff error. Defaults to `eps()` plus [`TINY`](@ref).

!!! warning
    Bounds can't be adaptively computed over nonstandard grids; passing a value for
    `adapt` greater than 0 when `m::Nonstandard` results in an error.

!!! note
    Calling [`FiniteDifferenceMethod`](@ref) objects is equivalent to passing them to `fdm`.

# Examples

```julia-repl
julia> fdm(central_fdm(5, 1), sin, 1; adapt=2)
0.5403023058681039

julia> fdm(central_fdm(2, 1), exp, 0, Val(true))
(FiniteDifferenceMethod:
  order of method:       2
  order of derivative:   1
  grid:                  [-1, 1]
  coefficients:          [-0.5, 0.5]
  roundoff error:        1.42e-14
  bounds on derivatives: 1.00e+02
  step size:             1.69e-08
  accuracy:              1.69e-06
, 1.0000000031817473)
```
"""
function fdm(
    m::M,
    f,
    x,
    ::Val{true};
    condition=DEFAULT_CONDITION,
    bound=_estimate_bound(f(x), condition),
    eps=(Base.eps(float(bound)) + TINY),
    adapt=m.history.adapt,
    max_step=0.1,
) where M<:FiniteDifferenceMethod
    if M <: Nonstandard && adapt > 0
        throw(ArgumentError("can't adaptively compute bounds over Nonstandard grids"))
    end
    eps > 0 || throw(ArgumentError("eps must be positive, got $eps"))
    bound > 0 || throw(ArgumentError("bound must be positive, got $bound"))
    0 <= adapt < 20 - m.p || throw(ArgumentError("can't perform $adapt adaptation steps"))

    # The below calculation can fail for `m.q == 0`, because then `C₂` may be zero. We
    # therefore handle this edge case here.
    m.q == 0 && return m, f(x)

    p = m.p
    q = m.q
    grid = m.grid
    coefs = m.coefs

    # Adaptively compute the bound on the function and derivative values, if applicable.
    if adapt > 0
        newm = (M.name.wrapper)(p + 1, p)
        dfdx = fdm(
            newm,
            f,
            x;
            condition=condition,
            eps=eps,
            bound=bound,
            max_step=max_step,
            adapt=(adapt - 1),
        )
        bound = _estimate_bound(dfdx, condition)
    end

    # Set the step size by minimising an upper bound on the error of the estimate.
    C₁ = eps * sum(abs, coefs)
    C₂ = bound * sum(n->abs(coefs[n] * grid[n]^p), eachindex(coefs)) / factorial(p)
    ĥ = min((q / (p - q) * C₁ / C₂)^(1 / p), max_step)

    # Estimate the accuracy of the method.
    accuracy = ĥ^(-q) * C₁ + ĥ^(p - q) * C₂

    # Estimate the value of the derivative.
    dfdx = sum(i->coefs[i] * f(x + ĥ * grid[i]), eachindex(grid)) / ĥ^q

    m.history.eps = eps
    m.history.bound = bound
    m.history.step = ĥ
    m.history.accuracy = accuracy

    return m, dfdx
end

function fdm(m::FiniteDifferenceMethod, f, x, ::Val{false}=Val(false); kwargs...)
    _, dfdx = fdm(m, f, x, Val(true); kwargs...)
    return dfdx
end

## Deprecations

# Used for keyword argument name deprecations
function _dep_kwarg(kwargs)
    for (old, new) in [(:ε, :eps), (:M, :bound)]
        haskey(kwargs, old) || continue
        val = kwargs[old]
        error(
            "keyword argument `", old, "` should now be passed as `", new, "` upon ",
            "application of the method. For example:\n    ",
            "central_fdm(5, 1)(f, x; $new=$val)\n",
            "not\n    ",
            "central_fdm(5, 1; $old=$val)(f, x)"
        )
    end
end

function fdm(
    grid::AbstractVector{<:Real},
    q::Int,
    ::Union{Val{true}, Val{false}}=Val(false);
    kwargs...,
)
    error("to use a custom grid, use `Nonstandard(grid, q)` and pass the result to `fdm`")
end

for fdmf in (:central_fdm, :backward_fdm, :forward_fdm)
    @eval function $fdmf(p::Int, q::Int, ::Union{Val{true}, Val{false}}; kwargs...)
        error(
            "the `Val` argument should now be passed directly to `fdm` after ",
            "constructing the method, not to the method constructor itself"
        )
    end
end


text:
export CRNN_Batch_TCS_With_Softmax
export DNN_Batch_TCS_With_Softmax
export RNN_Batch_TCS_With_Softmax

"""
    DNN_Batch_TCS_With_Softmax(x::Variable{Array{T}},
                               seqlabels::Vector,
                               inputlens;
                               background::Int=1,
                               foreground::Int=2,
                               weight=1.0) where T

a batch of concatenated input sequence is processed by neural networks into `x`

# Inputs
`x`         : 2-D Variable, a batch of concatenated input sequence\n
`seqlabels` : a batch of sequential labels, like [[i,j,k],[x,y],...]\n
`inputlens` : records each input sequence's length, like [20,17,...]\n
`weight`    : weight for TCS loss

# Structure
    ┌───┐
    │ │ │
    │ W ├──►─┐
    │ │ │    │
    └───┘    │
    ┌───┐    │    ┌───┐          ┌───┐
    │ │ │  ┌─┴─┐  │ │ │ softmax  │ │ │   ┌───────┐
    │ Z ├─►│ × ├─►│ X ├─────────►│ P ├──►│TCSLOSS│◄── (seqLabel)
    │ │ │  └───┘  │ │ │          │ │ │   └───┬───┘
    └───┘         └─┬─┘          └─┬─┘       │
                    │              │+        ▼
                  ┌─┴─┐            ▼       ┌─┴─┐
                  │ │ │          ┌─┴─┐ -   │ │ │
                  │ δ │◄─────────┤ - │──◄──┤ r │
                  │ │ │          └───┘     │ │ │
                  └───┘                    └───┘
"""
function DNN_Batch_TCS_With_Softmax(x::Variable{Array{T}},
                                    seqlabels::Vector,
                                    inputlens;
                                    background::Int=1,
                                    foreground::Int=2,
                                    weight=1.0) where T
    batchsize = length(seqlabels)
    loglikely = zeros(T, batchsize)
    I, F = indexbounds(inputlens)
    p = softmax(ᵛ(x); dims=1)
    r = zero(ᵛ(x))

    Threads.@threads for b = 1:batchsize
        span = I[b]:F[b]
        r[:,span], loglikely[b] = TCS(p[:,span], seqlabels[b], background=background, foreground=foreground)
        loglikely[b] /= length(seqlabels[b]) * 3 + 1
    end

    if x.backprop
        function DNN_Batch_TCS_With_Softmax_Backward()
            if need2computeδ!(x)
                if weight==1.0
                    δ(x) .+=  p - r
                else
                    δ(x) .+= (p - r) .* weight
                end
            end
        end
        push!(graph.backward, DNN_Batch_TCS_With_Softmax_Backward)
    end
    return sum(loglikely)/batchsize
end

"""
    RNN_Batch_TCS_With_Softmax(x::Variable{Array{T}},
                               seqlabels::Vector,
                               inputlens;
                               background::Int=1,
                               foreground::Int=2,
                               weight=1.0) where T

a batch of padded input sequence is processed by neural networks into `x`

# Inputs
`x`         : 3-D Variable with shape (featdims,timesteps,batchsize), a batch of padded input sequence\n
`seqlabels` : a batch of sequential labels, like [[i,j,k],[x,y],...]\n
`inputlens` : records each input sequence's length, like [20,17,...]\n
`weight`    : weight for TCS loss

# Structure
    ┌───┐
    │ │ │
    │ W ├──►─┐
    │ │ │    │
    └───┘    │
    ┌───┐    │    ┌───┐          ┌───┐
    │ │ │  ┌─┴─┐  │ │ │ softmax  │ │ │   ┌───────┐
    │ Z ├─►│ × ├─►│ X ├─────────►│ P ├──►│TCSLOSS│◄── (seqLabel)
    │ │ │  └───┘  │ │ │          │ │ │   └───┬───┘
    └───┘         └─┬─┘          └─┬─┘       │
                    │              │+        ▼
                  ┌─┴─┐            ▼       ┌─┴─┐
                  │ │ │          ┌─┴─┐ -   │ │ │
                  │ δ │◄─────────┤ - │──◄──┤ r │
                  │ │ │          └───┘     │ │ │
                  └───┘                    └───┘
"""
function RNN_Batch_TCS_With_Softmax(x::Variable{Array{T}},
                                    seqlabels::Vector,
                                    inputlens;
                                    background::Int=1,
                                    foreground::Int=2,
                                    weight=1.0) where T
    batchsize = length(seqlabels)
    loglikely = zeros(T, batchsize)
    p = zero(ᵛ(x))
    r = zero(ᵛ(x))

    Threads.@threads for b = 1:batchsize
        Tᵇ = inputlens[b]
        Lᵇ = length(seqlabels[b])
        p[:,1:Tᵇ,b] = softmax(x.value[:,1:Tᵇ,b]; dims=1)
        r[:,1:Tᵇ,b], loglikely[b] = TCS(p[:,1:Tᵇ,b], seqlabels[b], background=background, foreground=foreground)
        loglikely[b] /= Lᵇ * 3 + 1
    end

    if x.backprop
        function RNN_Batch_TCS_With_Softmax_Backward()
            if need2computeδ!(x)
                if weight==1.0
                    δ(x) .+=  p - r
                else
                    δ(x) .+= (p - r) .* weight
                end
            end
        end
        push!(graph.backward, RNN_Batch_TCS_With_Softmax_Backward)
    end
    return sum(loglikely)/batchsize
end

"""
    CRNN_Batch_TCS_With_Softmax(x::Variable{Array{T}},
                                seqlabels::Vector;
                                background::Int=1,
                                foreground::Int=2,
                                weight=1.0) where T

a batch of padded input sequence is processed by neural networks into `x`

# Main Inputs
`x`            : 3-D Variable with shape (featdims,timesteps,batchsize), resulted by a batch of padded input sequence\n
`seqlabels`    : a batch of sequential labels, like [[i,j,k],[x,y],...]\n
`weight`       : weight for TCS loss

# Structure
    ┌───┐
    │ │ │
    │ W ├──►─┐
    │ │ │    │
    └───┘    │
    ┌───┐    │    ┌───┐          ┌───┐
    │ │ │  ┌─┴─┐  │ │ │ softmax  │ │ │   ┌───────┐
    │ Z ├─►│ × ├─►│ X ├─────────►│ P ├──►│TCSLOSS│◄── (seqLabel)
    │ │ │  └───┘  │ │ │          │ │ │   └───┬───┘
    └───┘         └─┬─┘          └─┬─┘       │
                    │              │+        ▼
                  ┌─┴─┐            ▼       ┌─┴─┐
                  │ │ │          ┌─┴─┐ -   │ │ │
                  │ δ │◄─────────┤ - │──◄──┤ r │
                  │ │ │          └───┘     │ │ │
                  └───┘                    └───┘
"""
function CRNN_Batch_TCS_With_Softmax(x::Variable{Array{T}},
                                     seqlabels::Vector;
                                     background::Int=1,
                                     foreground::Int=2,
                                     weight=1.0) where T
    featdims, timesteps, batchsize = size(x)
    loglikely = zeros(T, batchsize)
    p = softmax(ᵛ(x); dims=1)
    r = zero(ᵛ(x))

    Threads.@threads for b = 1:batchsize
        r[:,:,b], loglikely[b] = TCS(p[:,:,b], seqlabels[b], background=background, foreground=foreground)
        loglikely[b] /= length(seqlabels[b]) * 3 + 1
    end

    if x.backprop
        function CRNN_Batch_TCS_With_Softmax_Backward()
            if need2computeδ!(x)
                if weight==1.0
                    δ(x) .+=  p - r
                else
                    δ(x) .+= (p - r) .* weight
                end
            end
        end
        push!(graph.backward, CRNN_Batch_TCS_With_Softmax_Backward)
    end
    return sum(loglikely)/batchsize
end

max_stars_repo_name: sonosole/Delta

max_issues_repo_name: sonosole/Delta

max_forks_repo_name: sonosole/Delta