text:
"""
    Initialize{GEBT} <: AbstractModel

Initial condition analysis model for the ``GEBT`` model.  This model may be used
to define a consistent set of initial conditions for use with the GEBT model.
State variables for this model are as defined by GXBeam.  Inputs correspond to
the external forces ``F_{x,i}, F_{y,i}, F_{z,i}, M_{x,i}, M_{y,i}, M_{z,i}`` or
displacements ``u_{x,i}, u_{y,i}, u_{z,i}, \\theta_{x,i}, \\theta_{y,i},
\\theta_{z,i}`` applied to each node, followed by the distributed loads
``f_{x,i}, f_{y,i}, f_{z,i}, m_{x,i}, m_{y,i}, m_{z,i}`` applied to each beam
element, followed by the linear and angular velocity of the system. Parameters
correspond to the location ``p_{x}, p_{y}, p_{z}`` of each node, followed by
each beam element's properties, followed by the initial displacement and velocity
of each beam element ``u_{x}, u_{y}, u_{z}, \dot{u}_x, \dot{u}_y, \dot{u}_z``.
Each beam element's properties are defined by a triad which defines the
orientation of the beam element ``e_{1,x}, e_{1,y}, e_{1,z}, e_{2,x}, e_{2,y},
e_{2,z}, e_{3,x}, e_{3,y}, e_{3,z}``, followed by the 21 independent entries of
the compliance matrix ``C_{11}, C_{12}, C_{13}, C_{14}, C_{15}, C_{16}, C_{22},
C_{23}, C_{24}, C_{25}, C_{26}, C_{33}, C_{34}, C_{35}, C_{36}, C_{44}, C_{45},
C_{46}, C_{55}, C_{56}, C_{66}``, followed by the beam element's inertial
properties ``\\mu, x_{m,2}, x_{m,3}, i_{22}, i_{33}, i_{23}``.
"""
struct Initialize{T}
    model::T
end

# --- Traits --- #

function number_of_states(model::Initialize)
    return number_of_states(model.model)
end

function number_of_inputs(model::Initialize)
    return number_of_inputs(model.model)
end

function number_of_parameters(model::Initialize{<:GEBT})
    return number_of_parameters(model.model) + 6*length(model.model.icol_elem)
end

inplaceness(::Type{Initialize{T}}) where T = inplaceness(T)

mass_matrix_type(::Type{Initialize{T}}) where T = mass_matrix_type(T)

state_jacobian_type(::Type{Initialize{T}}) where T = state_jacobian_type(T)

input_jacobian_type(::Type{Initialize{T}}) where T = input_jacobian_type(T)

# --- Methods --- #

function get_rates!(Mdq, model::Initialize{<:GEBT}, q, r, p, t)

    # number of points and elements
    np = length(model.model.icol_point)
    ne = length(model.model.icol_elem)

    # index of beginning and end of each beam element
    start = model.model.start
    stop = model.model.stop

    # scaling parameters
    force_scaling = model.model.force_scaling
    mass_scaling = model.model.mass_scaling

    # model pointers
    irow_point = model.model.irow_point
    irow_elem = model.model.irow_elem
    irow_elem1 = model.model.irow_elem1
    irow_elem2 = model.model.irow_elem2
    icol_point = model.model.icol_point
    icol_elem = model.model.icol_elem

    # location and degree of freedom of displacement constraints
    displacement = model.model.displacement

    # construct assembly from parameters
    assembly = gxbeam_assembly(p, np, ne, start, stop)

    # extract newly created beam elements
    elements = assembly.elements

    # construct point and element loads from inputs
    prescribed, distributed = gxbeam_loads(r, np, ne, displacement, elements)

    # origin, linear velocity, and angular velocity
    x0 = @SVector zeros(3)
    v0 = @SVector [r[end-5], r[end-4], r[end-3]]
    ω0 = @SVector [r[end-2], r[end-1], r[end]]

    # return mass matrix multiplied state rates
    return gxbeam_initialize_rates!(Mdq, q, assembly, prescribed, distributed,
        force_scaling, mass_scaling, irow_point, irow_elem,
        irow_elem1, irow_elem2, icol_point, icol_elem,
        x0, v0, ω0)
end

function get_mass_matrix!(M, model::Initialize, q, r, p, t)

    # number of points and elements
    np = length(model.icol_point)
    ne = length(model.icol_elem)

    # index of beginning and end of each beam element
    start = model.start
    stop = model.stop

    # scaling parameters
    force_scaling = model.force_scaling
    mass_scaling = model.mass_scaling

    # construct assembly from parameters
    assembly = gxbeam_assembly(p, np, ne, start, stop)

    # extract model pointers
    irow_point = model.irow_point
    irow_elem = model.irow_elem
    irow_elem1 = model.irow_elem1
    irow_elem2 = model.irow_elem2
    icol_point = model.icol_point
    icol_elem = model.icol_elem

    # return mass matrix
    return gxbeam_mass_matrix!(M, q, assembly, force_scaling, mass_scaling,
        irow_point, irow_elem, irow_elem1, irow_elem2, icol_point, icol_elem)
end

# --- Performance Overloads --- #

function get_state_jacobian!(J, model::Initialize, q, r, p, t)

    # number of points and elements
    np = length(model.icol_point)
    ne = length(model.icol_elem)

    # index of beginning and end of each beam element
    start = model.start
    stop = model.stop

    # scaling parameters
    force_scaling = model.force_scaling
    mass_scaling = model.mass_scaling

    # model pointers
    irow_point = model.irow_point
    irow_elem = model.irow_elem
    irow_elem1 = model.irow_elem1
    irow_elem2 = model.irow_elem2
    icol_point = model.icol_point
    icol_elem = model.icol_elem

    # location and degree of freedom of displacement constraints
    displacement = model.displacement

    # construct assembly from parameters
    assembly = gxbeam_assembly(p, np, ne, start, stop)

    # extract newly created beam elements
    elements = assembly.elements

    # construct point and element loads from inputs
    prescribed, distributed = gxbeam_loads(r, np, ne, displacement, elements)

    # origin, linear velocity, and angular velocity
    x0 = @SVector zeros(3)
    v0 = @SVector [r[end-5], r[end-4], r[end-3]]
    ω0 = @SVector [r[end-2], r[end-1], r[end]]

    # return jacobian of right hand side with respect to state variables
    return gxbeam_state_jacobian!(J, q, assembly, prescribed,
        distributed, force_scaling, mass_scaling, irow_point, irow_elem,
        irow_elem1, irow_elem2, icol_point, icol_elem, x0, v0, ω0)
end

# TODO: Update this function to include prescribed conditions and body velocities
# function get_input_jacobian(Jr, model::Initialize, q, r, p, t)
#
#     # extract system, assembly, and distributed loads
#     system = model.system
#     assembly = model.assembly
#     distributed = model.distributed
#
#     # extract start and stop of each beam element
#     elements = assembly.elements
#     start = assembly.start
#     stop = assembly.stop
#
#     # extract system constants and pointers
#     force_scaling = system.force_scaling
#     mass_scaling = system.mass_scaling
#     irow_point = system.irow_point
#
#     # get elements to which distributed loads are applied
#     element_indices = keys(distributed)
#
#     # return input jacobian linear map
#     return gxbeam_input_jacobian(q, r, elements, start, stop, force_scaling,
#         mass_scaling, irow_point, element_indices)
# end

# --- Unit Testing Methods --- #

function get_lhs(model::Initialize, dq, q, r, p, t)

    # number of points and elements
    np = length(model.icol_point)
    ne = length(model.icol_elem)

    # index of beginning and end of each beam element
    start = model.start
    stop = model.stop

    # scaling parameters
    force_scaling = model.force_scaling
    mass_scaling = model.mass_scaling

    # model pointers
    irow_point = model.irow_point
    irow_elem = model.irow_elem
    irow_elem1 = model.irow_elem1
    irow_elem2 = model.irow_elem2
    icol_point = model.icol_point
    icol_elem = model.icol_elem

    # location and degree of freedom of displacement constraints
    displacement = model.displacement

    # construct assembly from parameters
    assembly = gxbeam_assembly(p, np, ne, start, stop)

    # extract newly created beam elements
    elements = assembly.elements

    # construct point and element loads from inputs
    prescribed, distributed = gxbeam_loads(r, np, ne, displacement, elements)

    # origin, linear velocity, and angular velocity
    x0 = @SVector zeros(3)
    v0 = @SVector [r[end-5], r[end-4], r[end-3]]
    ω0 = @SVector [r[end-2], r[end-1], r[end]]

    return gxbeam_lhs(q, dq, assembly, prescribed, distributed, force_scaling,
        mass_scaling, irow_point, irow_elem, irow_elem1, irow_elem2, icol_point,
        icol_elem, x0, v0, ω0)
end

# --- Convenience Functions --- #

function set_states!(x, model::Initialize;
    u_e = nothing, theta_e = nothing, F_e = nothing, M_e = nothing,
    P_e = nothing, H_e = nothing, u_p = nothing, theta_p = nothing,
    F_p = nothing, M_p = nothing)

    icol_elem = model.icol_elem
    icol_point = model.icol_point

    displacement = model.displacement

    force_scaling = model.force_scaling
    mass_scaling = model.mass_scaling

    nelem = length(icol_elem)
    npoint = length(icol_point)

    for ielem = 1:nelem

        icol = icol_elem[ielem]

        if !isnothing(u_e)
            GXBeam.set_element_deflections!(x, icol, u_e[ielem])
        end

        if !isnothing(theta_e)
            GXBeam.set_element_rotations!(x, icol, theta_e[ielem])
        end

        if !isnothing(F_e)
            GXBeam.set_element_forces!(x, icol, F_e[ielem], force_scaling)
        end

        if !isnothing(M_e)
            GXBeam.set_element_moments!(x, icol, M_e[ielem], force_scaling)
        end

        if !isnothing(P_e)
            GXBeam.set_element_linear_momenta!(x, icol, P_e[ielem], mass_scaling)
        end

        if !isnothing(H_e)
            GXBeam.set_element_angular_momenta!(x, icol, H_e[ielem], mass_scaling)
        end
    end

    for ipoint = 1:npoint

        icol = icol_point[ipoint]

        prescribed_forces = SVector{6}(view(displacement, :, ipoint)) .== false

        if !isnothing(u_p)
            GXBeam.set_point_deflections!(x, icol, u_p[ipoint], prescribed_forces)
        end

        if !isnothing(theta_p)
            GXBeam.set_point_rotations!(x, icol, theta_p[ipoint], prescribed_forces)
        end

        if !isnothing(F_p)
            GXBeam.set_point_forces!(x, icol, F_p[ipoint], prescribed_forces, force_scaling)
        end

        if !isnothing(M_p)
            GXBeam.set_point_moments!(x, icol, M_p[ipoint], prescribed_forces, force_scaling)
        end
    end

    return x
end

function set_inputs!(y, model::Initialize; point_conditions, element_loads, V, Omega)

    np = length(model.icol_point)
    ne = length(model.icol_elem)

    for ip = 1:np
        y[6*(ip-1) + 1] = point_conditions[1,ip]
        y[6*(ip-1) + 2] = point_conditions[2,ip]
        y[6*(ip-1) + 3] = point_conditions[3,ip]
        y[6*(ip-1) + 4] = point_conditions[4,ip]
        y[6*(ip-1) + 5] = point_conditions[5,ip]
        y[6*(ip-1) + 6] = point_conditions[6,ip]
    end

    for ie = 1:ne
        y[6*np + 6*(ie-1) + 1] = element_loads[1,ie]
        y[6*np + 6*(ie-1) + 2] = element_loads[2,ie]
        y[6*np + 6*(ie-1) + 3] = element_loads[3,ie]
        y[6*np + 6*(ie-1) + 4] = element_loads[4,ie]
        y[6*np + 6*(ie-1) + 5] = element_loads[5,ie]
        y[6*np + 6*(ie-1) + 6] = element_loads[6,ie]
    end

    y[6*np + 6*ne + 1] = V[1]
    y[6*np + 6*ne + 2] = V[2]
    y[6*np + 6*ne + 3] = V[3]

    y[6*np + 6*ne + 4] = Omega[1]
    y[6*np + 6*ne + 5] = Omega[2]
    y[6*np + 6*ne + 6] = Omega[3]

    return y
end

function set_parameters!(p, model::Initialize; assembly)

    np = length(model.icol_point)
    ne = length(model.icol_elem)

    for ip = 1:np
        p[3*(ip-1) + 1] = assembly.points[ip][1]
        p[3*(ip-1) + 2] = assembly.points[ip][2]
        p[3*(ip-1) + 3] = assembly.points[ip][3]
    end

    for ie = 1:ne
        element = assembly.elements[ie]
        minv11, minv12, minv22 = element.minv11, element.minv12, element.minv22
        Msub = inv((@SMatrix [
            minv11[1,1] minv12[1,2] minv12[1,3];
            minv12[1,2] minv22[2,2] minv22[2,3];
            minv12[1,3] minv22[3,2] minv22[3,3];
            ]))
        p[3*np + 36*(ie-1) + 1] = element.Cab[1,1]
        p[3*np + 36*(ie-1) + 2] = element.Cab[2,1]
        p[3*np + 36*(ie-1) + 3] = element.Cab[3,1]
        p[3*np + 36*(ie-1) + 4] = element.Cab[1,2]
        p[3*np + 36*(ie-1) + 5] = element.Cab[2,2]
        p[3*np + 36*(ie-1) + 6] = element.Cab[3,2]
        p[3*np + 36*(ie-1) + 7] = element.Cab[1,3]
        p[3*np + 36*(ie-1) + 8] = element.Cab[2,3]
        p[3*np + 36*(ie-1) + 9] = element.Cab[3,3]
        p[3*np + 36*(ie-1) + 10] = element.C11[1,1]
        p[3*np + 36*(ie-1) + 11] = element.C11[1,2]
        p[3*np + 36*(ie-1) + 12] = element.C11[1,3]
        p[3*np + 36*(ie-1) + 13] = element.C12[1,1]
        p[3*np + 36*(ie-1) + 14] = element.C12[1,2]
        p[3*np + 36*(ie-1) + 15] = element.C12[1,3]
        p[3*np + 36*(ie-1) + 16] = element.C11[2,2]
        p[3*np + 36*(ie-1) + 17] = element.C11[2,3]
        p[3*np + 36*(ie-1) + 18] = element.C12[2,1]
        p[3*np + 36*(ie-1) + 19] = element.C12[2,2]
        p[3*np + 36*(ie-1) + 20] = element.C12[2,3]
        p[3*np + 36*(ie-1) + 21] = element.C11[3,3]
        p[3*np + 36*(ie-1) + 22] = element.C12[3,1]
        p[3*np + 36*(ie-1) + 23] = element.C12[3,2]
        p[3*np + 36*(ie-1) + 24] = element.C12[3,3]
        p[3*np + 36*(ie-1) + 25] = element.C22[1,1]
        p[3*np + 36*(ie-1) + 26] = element.C22[1,2]
        p[3*np + 36*(ie-1) + 27] = element.C22[1,3]
        p[3*np + 36*(ie-1) + 28] = element.C22[2,2]
        p[3*np + 36*(ie-1) + 29] = element.C22[2,3]
        p[3*np + 36*(ie-1) + 30] = element.C22[3,3]
        p[3*np + 36*(ie-1) + 31] =  Msub[1,1] # μ
        p[3*np + 36*(ie-1) + 32] = -Msub[1,3]/Msub[1,1] # xm2
        p[3*np + 36*(ie-1) + 33] =  Msub[1,2]/Msub[1,1] # xm3
        p[3*np + 36*(ie-1) + 34] =  Msub[2,2] # i22
        p[3*np + 36*(ie-1) + 35] =  Msub[3,3] # i33
        p[3*np + 36*(ie-1) + 36] = -Msub[2,3] # i23
    end

    return p
end

function separate_states(model::Initialize, x;
    point_conditions = fill(NaN, 6, length(model.icol_point)))

    TF = eltype(x)
    np = length(model.icol_point)
    ne = length(model.icol_elem)

    u_p = Vector{SVector{3,TF}}(undef, np)
    theta_p = Vector{SVector{3,TF}}(undef, np)
    F_p = Vector{SVector{3,TF}}(undef, np)
    M_p = Vector{SVector{3,TF}}(undef, np)

    u_e = Vector{SVector{3,TF}}(undef, ne)
    theta_e = Vector{SVector{3,TF}}(undef, ne)
    F_e = Vector{SVector{3,TF}}(undef, ne)
    M_e = Vector{SVector{3,TF}}(undef, ne)
    P_e = Vector{SVector{3,TF}}(undef, ne)
    H_e = Vector{SVector{3,TF}}(undef, ne)

    for ip = 1:np

        icol = model.icol_point[ip]

        prescribed_forces = SVector{6}(view(model.displacement, :, ip)) .== false

        # get the displacement and rotations of the point
        u_p[ip] = SVector(ifelse(prescribed_forces[1], x[icol  ], point_conditions[1,ip]),
                    ifelse(prescribed_forces[2], x[icol+1], point_conditions[2,ip]),
                    ifelse(prescribed_forces[3], x[icol+2], point_conditions[3,ip]))
        theta_p[ip] = SVector(ifelse(prescribed_forces[4], x[icol+3], point_conditions[4,ip]),
                    ifelse(prescribed_forces[5], x[icol+4], point_conditions[5,ip]),
                    ifelse(prescribed_forces[6], x[icol+5], point_conditions[6,ip]))

        # overwrite external forces/moments with solved for forces/moments
        F_p[ip] = SVector(ifelse(prescribed_forces[1], point_conditions[1,ip], x[icol  ] * model.force_scaling),
                    ifelse(prescribed_forces[2], point_conditions[2,ip], x[icol+1] * model.force_scaling),
                    ifelse(prescribed_forces[3], point_conditions[3,ip], x[icol+2] * model.force_scaling))
        M_p[ip] = SVector(ifelse(prescribed_forces[4], point_conditions[4,ip], x[icol+3] * model.force_scaling),
                    ifelse(prescribed_forces[5], point_conditions[5,ip], x[icol+4] * model.force_scaling),
                    ifelse(prescribed_forces[6], point_conditions[6,ip], x[icol+5] * model.force_scaling))

        # convert rotation parameter to Wiener-Milenkovic parameters
        scaling = GXBeam.rotation_parameter_scaling(theta_p[ip])
        theta_p[ip] *= scaling

    end

    for ie = 1:ne

        icol = model.icol_elem[ie]

        u_e[ie] = SVector(x[icol], x[icol+1], x[icol+2])
        theta_e[ie] = SVector(x[icol+3], x[icol+4], x[icol+5])
        F_e[ie] = SVector(x[icol+6], x[icol+7], x[icol+8]) .* model.force_scaling
        M_e[ie] = SVector(x[icol+9], x[icol+10], x[icol+11]) .* model.force_scaling
        P_e[ie] = SVector(x[icol+12], x[icol+13], x[icol+14]) .* model.mass_scaling
        H_e[ie] = SVector(x[icol+15], x[icol+16], x[icol+17]) .* model.mass_scaling

        # convert rotation parameter to Wiener-Milenkovic parameters
        scaling = GXBeam.rotation_parameter_scaling(theta_e[ie])
        theta_e[ie] *= scaling

    end

    return (u_p = u_p, theta_p = theta_p, F_p = F_p, M_p = M_p,
        u_e = u_e, theta_e = theta_e, F_e = F_e, M_e = M_e, P_e = P_e, H_e = H_e)
end

function separate_inputs(model::Initialize, y)

    np = length(model.icol_point)
    ne = length(model.icol_elem)

    point_conditions = reshape(view(y, 1 : 6*np), 6, np)

    element_loads = reshape(view(y, 6*np + 1 : 6*np + 6*ne), 6, ne)

    V = view(y, 6*np + 6*ne + 1 : 6*np + 6*ne + 3)

    Omega = view(y, 6*np + 6*ne + 4 : 6*np + 6*ne + 6)

    return (point_conditions = point_conditions, element_loads = element_loads,
        V = V, Omega = Omega)
end

function separate_parameters(model::Initialize, p)
    np = length(model.icol_point)
    ne = length(model.icol_elem)
    start = model.start
    stop = model.stop
    return (assembly = gxbeam_assembly(p, np, ne, start, stop),)
end

# --- Internal --- #

function gxbeam_assembly(p, np, ne, start, stop)
    points = [gxbeam_point(view(p, 3*(ip-1) + 1 : 3*(ip-1) + 3)) for ip = 1:np]
    elements = [gxbeam_element(view(p, 3*np + 36*(ie-1) + 1 : 3*np + 36*(ie-1) + 36),
        points, start[ie], stop[ie]) for ie = 1:ne]
    return GXBeam.Assembly(points, start, stop, elements)
end

gxbeam_point(p) = SVector(p[1], p[2], p[3])

function gxbeam_element(p, points, start, stop)
    # separate element parameters
    e1x, e1y, e1z, e2x, e2y, e2z, e3x, e3y, e3z, C11, C12, C13, C14, C15, C16,
        C22, C23, C24, C25, C26, C33, C34, C35, C36, C44, C45, C46, C55, C56,
        C66, μ, xm2, xm3, i22, i33, i23 = p
    # element length
    ΔL = norm(points[stop] - points[start])
    # element location
    x = (points[start] + points[stop])/2
    # element compliance matrix
    C = @SMatrix [
        C11 C12 C13 C14 C15 C16;
        C12 C22 C23 C24 C25 C26;
        C13 C23 C33 C34 C35 C36;
        C14 C24 C34 C44 C45 C46;
        C15 C25 C35 C45 C55 C56;
        C16 C26 C36 C46 C56 C66
        ]
    # element inverse mass matrix
    minv = inv(@SMatrix [
            μ       0     0       0 μ*xm3 -μ*xm2;
            0       μ     0  -μ*xm3     0      0;
            0       0     μ   μ*xm2     0      0;
            0  -μ*xm3 μ*xm2 i22+i33     0      0;
         μ*xm3      0     0       0   i22   -i23;
        -μ*xm2      0     0       0  -i23    i33;
    ])
    # element triad
    Cab = @SMatrix [
        e1x e2x e3x;
        e1y e2y e3y;
        e1z e2z e3z;
    ]

    return GXBeam.Element(ΔL, x, C, minv, Cab)
end

function gxbeam_loads(r, np, ne, d, elements)
    prescribed = Dict(ip => gxbeam_point_load(view(r,
        6*(ip-1) + 1 : 6*(ip-1) + 6), view(d, :, ip)) for ip = 1:np)
    distributed = Dict(ie => gxbeam_distributed_load(view(r,
        6*np + 6*(ie-1) + 1 : 6*np + 6*(ie-1) + 6), elements[ie].L) for ie = 1:ne)
    return prescribed, distributed
end

function gxbeam_point_load(r, d)
    displacement = SVector{6}(d)
    force = displacement .== false
    value = SVector{6}(r)
    follower = @SVector zeros(6)
    return PrescribedConditions(force, value, follower)
end

function gxbeam_distributed_load(r, ΔL)
    f = SVector(r[1], r[2], r[3])
    m = SVector(r[4], r[5], r[6])
    f1 = f2 = ΔL*f/2
    m1 = m2 = ΔL*m/2
    f1_follower = f2_follower = @SVector zeros(3)
    m1_follower = m2_follower = @SVector zeros(3)
    return DistributedLoads(f1, f2, m1, m2, f1_follower, f2_follower, m1_follower,
        m2_follower)
end

function gxbeam_lhs(q, dq, assembly, prescribed, distributed,
    force_scaling, mass_scaling, irow_point, irow_elem, irow_elem1, irow_elem2,
    icol_point, icol_elem, x0, v0, ω0)

    steady_residual = similar(dq)
    dynamic_residual = similar(dq)

    GXBeam.steady_state_system_residual!(steady_residual, q, assembly, prescribed,
        distributed, force_scaling, mass_scaling, irow_point, irow_elem, irow_elem1,
        irow_elem2, icol_point, icol_elem, x0, v0, ω0)

    GXBeam.dynamic_system_residual!(dynamic_residual, q, dq, assembly, prescribed,
        distributed, force_scaling, mass_scaling, irow_point, irow_elem, irow_elem1,
        irow_elem2, icol_point, icol_elem, x0, v0, ω0)

    j! = (J, x) -> system_jacobian!(J, x, assembly, pcond, dload, force_scaling,
        mass_scaling, irow_point, irow_elem, irow_elem1, irow_elem2, icol_point, icol_elem,
        x0, v0, ω0, u0, theta0, udot0, thetadot0)

    return steady_residual - dynamic_residual
end

function gxbeam_rhs!(out, u, assembly, prescribed,
    distributed, force_scaling, mass_scaling, irow_point, irow_elem,
    irow_elem1, irow_elem2, icol_point, icol_elem, x0, v0, ω0)

    GXBeam.system_residual!(F, x, assembly, pcond, dload, force_scaling,
        mass_scaling, irow_point, irow_elem, irow_elem1, irow_elem2, icol_point, icol_elem,
        x0, v0, ω0, u0, theta0, udot0, thetadot0)

    return out
end

function gxbeam_mass_matrix!(M, u, assembly, force_scaling, mass_scaling, irow_point,
    irow_elem, irow_elem1, irow_elem2, icol_point, icol_elem)

    M .= 0

    # mass matrix is GXBeam mass matrix, moved to LHS
    GXBeam.system_mass_matrix!(M, u, assembly, force_scaling, mass_scaling, irow_point,
        irow_elem, irow_elem1, irow_elem2, icol_point, icol_elem)

    M .*= -1

    return M
end

function gxbeam_rates!(Mdu, u, assembly, prescribed,
    distributed, force_scaling, mass_scaling, irow_point, irow_elem,
    irow_elem1, irow_elem2, icol_point, icol_elem, x0, v0, ω0)

    gxbeam_rhs!(Mdu, u, assembly, prescribed,
        distributed, force_scaling, mass_scaling, irow_point, irow_elem,
        irow_elem1, irow_elem2, icol_point, icol_elem, x0, v0, ω0)

    return Mdu
end

function gxbeam_state_jacobian!(K, u, assembly, prescribed,
    distributed, force_scaling, mass_scaling, irow_point, irow_elem,
    irow_elem1, irow_elem2, icol_point, icol_elem, x0, v0, ω0)

    K .= 0

    # jacobian of mass matrix multiplied rates is equal to steady state jacobian
    GXBeam.steady_state_system_jacobian!(K, u, assembly, prescribed,
        distributed, force_scaling, mass_scaling, irow_point, irow_elem, irow_elem1,
        irow_elem2, icol_point, icol_elem, x0, v0, ω0)

    return K
end

# function gxbeam_input_jacobian(q, r, elements, start, stop, force_scaling, mass_scaling,
#     irow_point, element_indices)
#
#     f! = (y, x) -> gxbeam_jac_vec!(y, x, r, elements, start, stop, force_scaling,
#         mass_scaling, irow_point, element_indices)
#     M = length(q)
#     N = length(r)
#
#     return LinearMap(f!, M, N; ismutating = true)
# end
#
# # input jacobian vector product
# function gxbeam_jac_vec!(y, x, r, elements, start, stop, force_scaling, mass_scaling,
#     irow_point, element_indices)
#
#     nelem = length(start)
#
#     # initialize output array
#     y .= 0
#
#     # current input array index
#     ir = 0
#
#     # distributed loads
#     for ielem = 1:nelem
#         if ielem in element_indices
#             # get point indices
#             irow_p1 = irow_point[start[ielem]]
#             irow_p2 = irow_point[stop[ielem]]
#             # initialize element load jacobians
#             f_r = Diagonal((@SVector ones(3)))
#             m_r = Diagonal((@SVector ones(3)))
#             # jacobians after integration
#             # (assume distributed loads are constant on each element)
#             f1_r = f2_r = elements[ielem].L*f_r/2
#             m1_r = m2_r = elements[ielem].L*m_r/2
#             # calculate jacobian vector product
#             xf = SVector(x[ir+1], x[ir+2], x[ir+3])
#             xm = SVector(x[ir+4], x[ir+5], x[ir+6])
#             y[irow_p1:irow_p1+2] .+= -f1_r * xf ./ force_scaling
#             y[irow_p1+3:irow_p1+5] .+= -m1_r * xm ./ force_scaling
#             y[irow_p2:irow_p2+2] .+= -f2_r * xf ./ force_scaling
#             y[irow_p2+3:irow_p2+5] .+= -m2_r * xm ./ force_scaling
#             # move to next set of distributed loads
#             ir += 6
#         end
#     end
#
#     return y
# end
#
# # update distributed loads using inputs
# function update_distributed(distributed, elements, r::AbstractVector{T}) where T
#     new_distributed = Dict{Int, DistributedLoads{T}}()
#     for (key, value) in distributed
#         new_distributed[key] = DistributedLoads{T}(value)
#     end
#     update_distributed!(new_distributed, elements, r)
#     return new_distributed
# end
#
# # update distributed loads using inputs
# function update_distributed!(distributed, elements, r)
#     distributed_points = keys(distributed)
#     largest_point = maximum(distributed_points)
#     ir = 0
#     for ielem in 1:largest_point
#         if ielem in distributed_points
#             # extract loads on this element
#             f = SVector(r[ir+1], r[ir+2], r[ir+3])
#             m = SVector(r[ir+4], r[ir+5], r[ir+6])
#             # get integrated distributed loads
#             f1 = f2 = elements[ielem].L*f/2
#             m1 = m2 = elements[ielem].L*m/2
#             # keep the same follower loads
#             f1_follower = distributed[ielem].f1_follower
#             f2_follower = distributed[ielem].f2_follower
#             m1_follower = distributed[ielem].m1_follower
#             m2_follower = distributed[ielem].m2_follower
#             # update distributed loads
#             distributed[ielem] = DistributedLoads(f1, f2, m1, m2, f1_follower,
#                 f2_follower, m1_follower, m2_follower)
#             ir += 6
#         end
#     end
#
#     return distributed
# end


text:
module TestSKL

using Random
using TSML
using TSMLextra
using PyCall
using Test

const IRIS = getiris()
const X = IRIS[:,1:4] |> DataFrame
const Y = IRIS[:,5] |> Vector

const XX = IRIS[:,2:4] |> DataFrame
const YY = IRIS[:,1] |> Vector

const classifiers = [
    "LinearSVC","QDA","MLPClassifier","BernoulliNB",
    "RandomForestClassifier",
    "NearestCentroid","SVC","LinearSVC","NuSVC","MLPClassifier",
    "SGDClassifier","KNeighborsClassifier",
    "DecisionTreeClassifier",
    "PassiveAggressiveClassifier","RidgeClassifier",
    "ExtraTreesClassifier","GradientBoostingClassifier",
    "BaggingClassifier","AdaBoostClassifier","GaussianNB","MultinomialNB",
    "ComplementNB","BernoulliNB"
    #"RidgeClassifierCV","LDA",
    #"GaussianProcessClassifier",
 ]

const regressors = [
    "SVR",
    "Ridge",
    "Lasso",
    "ElasticNet",
    "Lars",
    "LassoLars",
    "OrthogonalMatchingPursuit",
    "SGDRegressor",
    "PassiveAggressiveRegressor",
    "KNeighborsRegressor",
    "RadiusNeighborsRegressor",
    "DecisionTreeRegressor",
    "RandomForestRegressor",
    "ExtraTreesRegressor",
    "GradientBoostingRegressor",
    "MLPRegressor",
    "AdaBoostRegressor"
    #"RidgeCV",
    #"BayesianRidge",
    #"ARDRegression",
    #"KernelRidge",
    #"GaussianProcessRegressor",
    #"IsotonicRegression",
]

function fit_test(learner::String,in::DataFrame,out::Vector)
	_learner=SKLearner(Dict(:learner=>learner))
	fit!(_learner,in,out)
	@test _learner.model != nothing
	return _learner
end

function fit_transform_reg(model::Learner,in::DataFrame,out::Vector)
    @test sum((transform!(model,in) .- out).^2)/length(out) < 10.0
end

@testset "scikit classifiers" begin
    Random.seed!(1)
    for cl in classifiers
	println(cl)
	fit_test(cl,X,Y)
    end
end

@testset "scikit regressors" begin
    Random.seed!(1)
    for rg in regressors
	println(rg)
	model=fit_test(rg,XX,YY)
	fit_transform_reg(model,XX,YY)
    end
end

end

ext: jl

max_stars_repo_name: ppalmes/TSMLextra

max_issues_repo_name: ppalmes/TSMLextra

max_forks_repo_name: ppalmes/TSMLextra


text:
struct ConnectionError <: Exception
  msg::String
end

function Base.showerror(io::IO, err::ConnectionError)
  print(io, "ConnectionError: ", err.msg)
end


text:
@testset "Features" begin

for l in 7:5:20
    for offset in 2π .* [0, 2, -1, 10, -50]
        phase_uw = collect(range(-2π; stop=2π, length=l))
        phase = rem2pi.(phase_uw, RoundNearest) .+ offset
        for mag in [nothing, collect(range(10; stop=1, length=l)), collect(range(1; stop=10, length=l))]
            unwrapped = unwrap(phase; mag=mag, correctglobal=true)
            @test unwrapped ≈ phase_uw
        end
    end
end

end


text:
# This file is a part of Julia. License is MIT: https://julialang.org/license

module IOTests

using Test
using Dates

const BASE_TEST_PATH = joinpath(Sys.BINDIR, "..", "share", "julia", "test")
include(joinpath(BASE_TEST_PATH, "testhelpers", "withlocales.jl"))

@testset "string/show representation of Date" begin
    @test string(DatesPlus.Date(1, 1, 1)) == "0001-01-01" # January 1st, 1 AD/CE
    @test sprint(show, DatesPlus.Date(1, 1, 1)) == "DatesPlus.Date(\"0001-01-01\")"
    @test string(DatesPlus.Date(0, 12, 31)) == "0000-12-31" # December 31, 1 BC/BCE
    @test DatesPlus.Date(1, 1, 1) - DatesPlus.Date(0, 12, 31) == DatesPlus.Day(1)
    @test DatesPlus.Date(DatesPlus.UTD(-306)) == DatesPlus.Date(0, 2, 29)
    @test string(DatesPlus.Date(0, 1, 1)) == "0000-01-01" # January 1st, 1 BC/BCE
    @test string(DatesPlus.Date(-1, 1, 1)) == "-0001-01-01" # January 1st, 2 BC/BCE
    @test string(DatesPlus.Date(-1000000, 1, 1)) == "-1000000-01-01"
    @test string(DatesPlus.Date(1000000, 1, 1)) == "1000000-01-01"
    @test string(DatesPlus.DateTime(2000, 1, 1, 0, 0, 0, 1)) == "2000-01-01T00:00:00.001"
    @test sprint(show, DatesPlus.DateTime(2000, 1, 1, 0, 0, 0, 1)) == "DatesPlus.DateTime(\"2000-01-01T00:00:00.001\")"
    @test string(DatesPlus.DateTime(2000, 1, 1, 0, 0, 0, 2)) == "2000-01-01T00:00:00.002"
    @test string(DatesPlus.DateTime(2000, 1, 1, 0, 0, 0, 500)) == "2000-01-01T00:00:00.500"
    @test string(DatesPlus.DateTime(2000, 1, 1, 0, 0, 0, 998)) == "2000-01-01T00:00:00.998"
    @test string(DatesPlus.DateTime(2000, 1, 1, 0, 0, 0, 999)) == "2000-01-01T00:00:00.999"
end

@testset "string/show representation of Time" begin
    tests = [
        DatesPlus.Time(0) => ("00:00:00", "DatesPlus.Time(0)"),
        DatesPlus.Time(0, 1) => ("00:01:00", "DatesPlus.Time(0, 1)"),
        DatesPlus.Time(0, 1, 2) => ("00:01:02", "DatesPlus.Time(0, 1, 2)"),
        DatesPlus.Time(0, 1, 2, 3) => ("00:01:02.003", "DatesPlus.Time(0, 1, 2, 3)"),
        DatesPlus.Time(0, 1, 2, 3, 4) => ("00:01:02.003004", "DatesPlus.Time(0, 1, 2, 3, 4)"),
        DatesPlus.Time(0, 1, 2, 3, 4, 5) => ("00:01:02.003004005", "DatesPlus.Time(0, 1, 2, 3, 4, 5)"),
        DatesPlus.Time(0, 0, 0, 0, 0, 1) => ("00:00:00.000000001", "DatesPlus.Time(0, 0, 0, 0, 0, 1)"),
        DatesPlus.Time(0, 0, 0, 1) => ("00:00:00.001", "DatesPlus.Time(0, 0, 0, 1)"),
    ]

    for (t, (printed, shown)) in tests
        @test sprint(print, t) == printed
        @test string(t) == printed
        @test sprint(show, t) == shown
        @test repr(t) == shown
    end
end

@testset "DateTime parsing" begin
    # Useful reference for different locales: http://library.princeton.edu/departments/tsd/katmandu/reference/months.html

    # Allow parsing of strings which are not representable as a TimeType
    str = "02/15/1996 25:00"
    df = DatesPlus.DateFormat("mm/dd/yyyy HH:MM")
    parsed = Any[
        DatesPlus.Month(2), DatesPlus.Day(15), DatesPlus.Year(1996), DatesPlus.Hour(25), DatesPlus.Minute(0)
    ]
    @test DatesPlus.parse_components(str, df) == parsed
    @test_throws ArgumentError DatesPlus.parse(DateTime, str, df)
end

@testset "DateFormat printing" begin
    @test sprint(show, DateFormat("yyyzzxmmdd\\MHH:MM:SS\\P")) == "dateformat\"yyyzzxmmdd\\MHH:MM:SSP\""
    @test sprint(show, DateFormat("yyy").tokens[1]) == "DatePart(yyy)"
    @test sprint(show, DateFormat("mmzzdd").tokens[2]) == "Delim(zz)"
    @test sprint(show, DateFormat("ddxmm").tokens[2]) == "Delim(x)"
    @test sprint(show, DateFormat("xxmmxx").tokens[2]) == "DatePart(mm)"
end

@testset "Common Parsing Patterns" begin
    #'1996-January-15'
    dt = DatesPlus.DateTime(1996, 1, 15)
    f = "yy-mm-dd"
    a = "96-01-15"
    @test DateTime(a, f) + DatesPlus.Year(1900) == dt
    @test DatesPlus.format(dt, f) == a
    a1 = "96-1-15"
    @test DatesPlus.DateTime(a1, f) + DatesPlus.Year(1900) == dt
    @test DatesPlus.format(dt, "yy-m-dd") == a1
    a2 = "96-1-1"
    @test DatesPlus.DateTime(a2, f) + DatesPlus.Year(1900) + DatesPlus.Day(14) == dt
    @test DatesPlus.format(dt - DatesPlus.Day(14), "yy-m-d") == a2
    a3 = "1996-1-15"
    @test DatesPlus.DateTime(a3, f) == dt
    @test DatesPlus.format(dt, "yyyy-m-d") == a3
    a4 = "1996-Jan-15"
    @test_throws ArgumentError DatesPlus.DateTime(a4, f) # Trying to use month name, but specified only "mm"

    f = "yy/uuu/dd"
    b = "96/Feb/15"
    @test DatesPlus.DateTime(b, f) + DatesPlus.Year(1900) == dt + DatesPlus.Month(1)
    @test DatesPlus.format(dt + DatesPlus.Month(1), f) == b
    b1 = "1996/Feb/15"
    @test DatesPlus.DateTime(b1, f) == dt + DatesPlus.Month(1)
    @test DatesPlus.format(dt + DatesPlus.Month(1), "yyyy/uuu/dd") == b1
    b2 = "96/Feb/1"
    @test DatesPlus.DateTime(b2, f) + DatesPlus.Year(1900) + DatesPlus.Day(14) == dt + DatesPlus.Month(1)
    @test DatesPlus.format(dt + DatesPlus.Month(1) - DatesPlus.Day(14), "yy/uuu/d") == b2
    # Here we've specified a text month name, but given a number
    b3 = "96/2/15"
    @test_throws ArgumentError DatesPlus.DateTime(b3, f)
    try
        DatesPlus.parse(DateTime, "2012/2/20T9:9:31.25i90", dateformat"yyyy/mm/ddTHH:MM:SS.s")
        @test false
    catch err
        @test isa(err, ArgumentError)
        @test err.msg == "Found extra characters at the end of date time string"
    end
    try
        DatesPlus.parse(DateTime, "2012/2/20T9:9:3i90", dateformat"yyyy/mm/ddTHH:MM:SS.s")
        @test false
    catch err
        @test isa(err, ArgumentError)
        @test err.msg == "Unable to parse date time. Expected directive Delim(.) at char 16"
    end

    f = "yy:dd:mm"
    c = "96:15:01"
    @test DatesPlus.DateTime(c, f) + DatesPlus.Year(1900) == dt
    @test DatesPlus.format(dt, f) == c
    c1 = "1996:15:01"
    @test DatesPlus.DateTime(c1, f) == dt
    @test DatesPlus.format(dt, "yyyy:dd:mm") == c1
    c2 = "96:15:1"
    @test DatesPlus.DateTime(c2, f) + DatesPlus.Year(1900) == dt
    @test DatesPlus.format(dt, "yy:dd:m") == c2
    c3 = "96:1:01"
    @test DatesPlus.DateTime(c3, f) + DatesPlus.Year(1900) + DatesPlus.Day(14) == dt
    @test DatesPlus.format(dt - DatesPlus.Day(14), "yy:m:dd") == c3
    c4 = "1996:15:01 # random comment"
    @test_throws ArgumentError DatesPlus.DateTime(c4, f)

    f = "yyyy, uuu, dd"
    d = "1996, Jan, 15"
    @test DatesPlus.DateTime(d, f) == dt
    @test DatesPlus.format(dt, f) == d
    d1 = "96, Jan, 15"
    @test DatesPlus.DateTime(d1, f) + DatesPlus.Year(1900) == dt
    @test DatesPlus.format(dt, "yy, uuu, dd") == d1
    d2 = "1996, Jan, 1"
    @test DatesPlus.DateTime(d2, f) + DatesPlus.Day(14) == dt
    @test DatesPlus.format(dt - DatesPlus.Day(14), "yyyy, uuu, d") == d2
    d3 = "1996, 2, 15"
    @test_throws ArgumentError DatesPlus.DateTime(d3, f)

    f = "yyyy.U.dd"
    e = "1996.January.15"
    @test DatesPlus.DateTime(e, f) == dt
    @test DatesPlus.format(dt, f) == e
    e1 = "96.January.15"
    @test DatesPlus.DateTime(e1, f) + DatesPlus.Year(1900) == dt
    @test DatesPlus.format(dt, "yy.U.dd") == e1

    fo = "yyyy m dd"
    f = "1996 1 15"
    @test DatesPlus.DateTime(f, fo) == dt
    @test DatesPlus.format(dt, fo) == f
    f1 = "1996 01 15"
    @test DatesPlus.DateTime(f1, fo) == dt
    @test DatesPlus.format(dt, "yyyy mm dd") == f1
    f2 = "1996 1 1"
    @test DatesPlus.DateTime(f2, fo) + DatesPlus.Day(14) == dt
    @test DatesPlus.format(dt - DatesPlus.Day(14), "yyyy m d") == f2

    j = "1996-01-15"
    f = "yyyy-mm-dd zzz"
    @test DatesPlus.DateTime(j, f) == dt
    @test DatesPlus.format(dt, f) == j * " zzz"
    k = "1996-01-15 10:00:00"
    f = "yyyy-mm-dd HH:MM:SS zzz"
    @test DatesPlus.DateTime(k, f) == dt + DatesPlus.Hour(10)
    @test DatesPlus.format(dt + DatesPlus.Hour(10), f) == k * " zzz"
    l = "1996-01-15 10:10:10.25"
    f = "yyyy-mm-dd HH:MM:SS.ss zzz"
    @test DatesPlus.DateTime(l, f) == dt + DatesPlus.Hour(10) + DatesPlus.Minute(10) + DatesPlus.Second(10) + DatesPlus.Millisecond(250)
    @test DatesPlus.format(dt + DatesPlus.Hour(10) + DatesPlus.Minute(10) + DatesPlus.Second(10) + DatesPlus.Millisecond(250), f) == l * " zzz"

    r = "1/15/1996" # Excel
    f = "m/dd/yyyy"
    @test DatesPlus.DateTime(r, f) == dt
    @test DatesPlus.format(dt, f) == r
    s = "19960115"
    f = "yyyymmdd"
    @test DatesPlus.DateTime(s, f) == dt
    @test DatesPlus.format(dt, f) == s
    v = "1996-01-15 10:00:00"
    f = "yyyy-mm-dd HH:MM:SS"
    @test DatesPlus.DateTime(v, f) == dt + DatesPlus.Hour(10)
    @test DatesPlus.format(dt + DatesPlus.Hour(10), f) == v
    w = "1996-01-15T10:00:00"
    f = "yyyy-mm-ddTHH:MM:SS zzz"
    @test DatesPlus.DateTime(w, f) == dt + DatesPlus.Hour(10)
    @test DatesPlus.format(dt + DatesPlus.Hour(10), f) == w * " zzz"

    f = "yyyy/m"
    y = "1996/1"
    @test DatesPlus.DateTime(y, f) == dt - DatesPlus.Day(14)
    @test DatesPlus.format(dt, f) == y
    y1 = "1996/1/15"
    @test_throws ArgumentError DatesPlus.DateTime(y1, f)
    y2 = "96/1"
    @test DatesPlus.DateTime(y2, f) + DatesPlus.Year(1900) == dt - DatesPlus.Day(14)
    @test DatesPlus.format(dt, "yy/m") == y2

    f = "yyyy"
    z = "1996"
    @test DatesPlus.DateTime(z, f) == dt - DatesPlus.Day(14)
    @test DatesPlus.format(dt, f) == z
    z1 = "1996-3"
    @test_throws ArgumentError DatesPlus.DateTime(z1, f)
    z2 = "1996-3-1"
    @test_throws ArgumentError DatesPlus.DateTime(z2, f)

    aa = "1/5/1996"
    f = "m/d/yyyy"
    @test DatesPlus.DateTime(aa, f) == dt - DatesPlus.Day(10)
    @test DatesPlus.format(dt - DatesPlus.Day(10), f) == aa
    bb = "5/1/1996"
    f = "d/m/yyyy"
    @test DatesPlus.DateTime(bb, f) == dt - DatesPlus.Day(10)
    @test DatesPlus.format(dt - DatesPlus.Day(10), f) == bb
    cc = "01151996"
    f = "mmddyyyy"
    @test DatesPlus.DateTime(cc, f) == dt
    @test DatesPlus.format(dt, f) == cc
    dd = "15011996"
    f = "ddmmyyyy"
    @test DatesPlus.DateTime(dd, f) == dt
    @test DatesPlus.format(dt, f) == dd
    ee = "01199615"
    f = "mmyyyydd"
    @test DatesPlus.DateTime(ee, f) == dt
    @test DatesPlus.format(dt, f) == ee
    ff = "1996-15-Jan"
    f = "yyyy-dd-uuu"
    @test DatesPlus.DateTime(ff, f) == dt
    @test DatesPlus.format(dt, f) == ff
    gg = "Jan-1996-15"
    f = "uuu-yyyy-dd"
    @test DatesPlus.DateTime(gg, f) == dt
    @test DatesPlus.format(dt, f) == gg
    hh = "1996#1#15"
    f = "yyyy#m#d"
    @test DatesPlus.DateTime(hh, f) == dt
    @test DatesPlus.format(dt, f) == hh

    f = "ymd"
    @test DatesPlus.Date("111", f) == DatesPlus.Date(1)
    @test DatesPlus.Date("1", f) == DatesPlus.Date(1)

    @test DatesPlus.DateTime("20140529 120000", "yyyymmdd HHMMSS") == DatesPlus.DateTime(2014, 5, 29, 12)

    @test DatesPlus.Date(string(DatesPlus.Date(dt))) == DatesPlus.Date(dt)
    @test DatesPlus.DateTime(string(dt)) == dt
end

@testset "prefix." begin
    s = "/1996/1/15"
    f = "/yyyy/m/d"
    dt = DatesPlus.DateTime(1996, 1, 15)
    @test DatesPlus.DateTime(s, f) == dt
    @test DatesPlus.format(dt, f) == s
    @test_throws ArgumentError DatesPlus.DateTime("1996/1/15", f)
end

@testset "French and Chinese" begin
    # from Jiahao
    @test DatesPlus.Date("2009年12月01日", "yyyy年mm月dd日") == DatesPlus.Date(2009, 12, 1)
    @test DatesPlus.format(DatesPlus.Date(2009, 12, 1), "yyyy年mm月dd日") == "2009年12月01日"
    @test DatesPlus.Date("2009-12-01", "yyyy-mm-dd") == DatesPlus.Date(2009, 12, 1)

    # French: from Milan
    f = "dd/mm/yyyy"
    f2 = "dd/mm/yy"
    @test DatesPlus.Date("28/05/2014", f) == DatesPlus.Date(2014, 5, 28)
    @test DatesPlus.Date("28/05/14", f2) + DatesPlus.Year(2000) == DatesPlus.Date(2014, 5, 28)
end

@testset "Customizing locale" begin
    DatesPlus.LOCALES["french"] = DatesPlus.DateLocale(
        ["janvier", "février", "mars", "avril", "mai", "juin",
         "juillet", "août", "septembre", "octobre", "novembre", "décembre"],
        ["janv", "févr", "mars", "avril", "mai", "juin",
         "juil", "août", "sept", "oct", "nov", "déc"],
        ["lundi", "mardi", "mercredi", "jeudi", "vendredi", "samedi", "dimanche"],
        [""],
    )

    f = "dd uuuuu yyyy"
    @test DatesPlus.Date("28 mai 2014", f; locale="french") == DatesPlus.Date(2014, 5, 28)
    @test DatesPlus.format(DatesPlus.Date(2014, 5, 28), f; locale="french") == "28 mai 2014"
    @test DatesPlus.Date("28 févr 2014", f; locale="french") == DatesPlus.Date(2014, 2, 28)
    @test DatesPlus.format(DatesPlus.Date(2014, 2, 28), f; locale="french") == "28 févr 2014"
    @test DatesPlus.Date("28 août 2014", f; locale="french") == DatesPlus.Date(2014, 8, 28)
    @test DatesPlus.format(DatesPlus.Date(2014, 8, 28), f; locale="french") == "28 août 2014"
    @test DatesPlus.Date("28 avril 2014", f; locale="french") == DatesPlus.Date(2014, 4, 28)
    @test DatesPlus.format(DatesPlus.Date(2014, 4, 28), f; locale="french") == "28 avril 2014"

    f = "dd u yyyy"
    @test DatesPlus.Date("28 avril 2014", f; locale="french") == DatesPlus.Date(2014, 4, 28)
    f = "dduuuuyyyy"
    # parses 3 and 4 character month names
    @test DatesPlus.Date("28mai2014", f; locale="french") == DatesPlus.Date(2014, 5, 28)
    @test DatesPlus.Date("28août2014", f; locale="french") == DatesPlus.Date(2014, 8, 28)
    # doesn't parse month name greater than 4 chars
    @test_throws ArgumentError DatesPlus.Date("28avril2014", f; locale="french")
end

@testset "year digits parsing" begin
    # From Tony Fong
    f = "dduuuyy"
    @test DatesPlus.Date("01Dec09", f) + DatesPlus.Year(2000) == DatesPlus.Date(2009, 12, 1)
    @test DatesPlus.format(DatesPlus.Date(2009, 12, 1), f) == "01Dec09"
    f = "dduuuyyyy"
    @test DatesPlus.Date("01Dec2009", f) == DatesPlus.Date(2009, 12, 1)
    @test DatesPlus.format(DatesPlus.Date(2009, 12, 1), f) == "01Dec2009"
    f = "duy"
    globex = ["f", "g", "h", "j", "k", "m", "n", "q", "u", "v", "x", "z"]
    locale = DatesPlus.DateLocale(globex, map(uppercase, globex), globex[1:7], globex[1:7])
    @test DatesPlus.Date("1F4", f; locale=locale) + DatesPlus.Year(2010) == DatesPlus.Date(2014, 1, 1)
    @test DatesPlus.format(DatesPlus.Date(2014, 1, 1), f; locale=locale) == "1F4"

    # From Matt Bauman
    f = "yyyy-mm-ddTHH:MM:SS"
    @test DatesPlus.DateTime("2014-05-28T16:46:04", f) == DatesPlus.DateTime(2014, 5, 28, 16, 46, 04)
end

@testset "Error handling" begin
    # Specified mm/dd, but date string has day/mm
    @test_throws ArgumentError DatesPlus.DateTime("18/05/2009", "mm/dd/yyyy")
    @test_throws ArgumentError DatesPlus.DateTime("18/05/2009 16", "mm/dd/yyyy hh")
    # Used "mm" for months AND minutes
    @test_throws ArgumentError DatesPlus.DateTime("18/05/2009 16:12", "mm/dd/yyyy hh:mm")
    # Date string has different delimiters than format string
    @test_throws ArgumentError DatesPlus.DateTime("18:05:2009", "mm/dd/yyyy")

    f = "y m d"
    @test DatesPlus.Date("1 1 1", f) == DatesPlus.Date(1)
    @test DatesPlus.Date("10000000000 1 1", f) == DatesPlus.Date(10000000000)
    @test_throws ArgumentError DatesPlus.Date("1 13 1", f)
    @test_throws ArgumentError DatesPlus.Date("1 1 32", f)
    @test_throws ArgumentError DatesPlus.Date(" 1 1 32", f)
    @test_throws ArgumentError DatesPlus.Date("# 1 1 32", f)
    @test DatesPlus.Date("1", f) == DatesPlus.Date(1)
    @test DatesPlus.Date("1 ", f) == DatesPlus.Date(1)
    @test DatesPlus.Date("1 2", f) == DatesPlus.Date(1, 2)
    # can't find space delimiter (finds '/'), so fails
    @test_throws ArgumentError DatesPlus.Date("2000/1", f)
end

@testset "formerly vectorized Date/DateTime/format methods" begin
    dr = ["2000-01-01", "2000-01-02", "2000-01-03", "2000-01-04", "2000-01-05",
          "2000-01-06", "2000-01-07", "2000-01-08", "2000-01-09", "2000-01-10"]
    dr2 = [DatesPlus.Date(2000) : DatesPlus.Day(1) : DatesPlus.Date(2000, 1, 10);]
    @test DatesPlus.Date.(dr) == dr2
    @test DatesPlus.Date.(dr, dateformat"yyyy-mm-dd") == dr2
    @test DatesPlus.DateTime.(dr) == DatesPlus.DateTime.(dr2)
    @test DatesPlus.DateTime.(dr, dateformat"yyyy-mm-dd") == DatesPlus.DateTime.(dr2)

    @test DatesPlus.format.(dr2, "yyyy-mm-dd") == dr

    @test typeof(DatesPlus.Date.(dr)) == Array{Date, 1}
end

@testset "Issue 13" begin
    t = DatesPlus.DateTime(1, 1, 1, 14, 51, 0, 118)
    @test DatesPlus.DateTime("[14:51:00.118]", "[HH:MM:SS.sss]") == t
    @test DatesPlus.DateTime("14:51:00.118", "HH:MM:SS.sss") == t
    @test DatesPlus.DateTime("[14:51:00.118?", "[HH:MM:SS.sss?") == t
    @test DatesPlus.DateTime("?14:51:00.118?", "?HH:MM:SS.sss?") == t
    @test DatesPlus.DateTime("x14:51:00.118", "xHH:MM:SS.sss") == t
    @test DatesPlus.DateTime("14:51:00.118]", "HH:MM:SS.sss]") == t
end

@testset "RFC1123Format" begin
    dt = DatesPlus.DateTime(2014, 8, 23, 17, 22, 15)
    @test DatesPlus.format(dt, DatesPlus.RFC1123Format) == "Sat, 23 Aug 2014 17:22:15"
    @test DatesPlus.DateTime(DatesPlus.format(dt, DatesPlus.RFC1123Format), DatesPlus.RFC1123Format) == dt
    @test DatesPlus.format(dt, "yyyy-mm-ddTHH:MM:SS E") == "2014-08-23T17:22:15 Saturday"
    @test DatesPlus.format(dt, "yyyy-mm-ddTHH:MM:SS e") == "2014-08-23T17:22:15 Sat"
    @test DatesPlus.format(dt, "yyyy-mm-dd E") == "2014-08-23 Saturday"
    @test DatesPlus.format(dt, "yyyy-mm-dd e") == "2014-08-23 Sat"
    @test DatesPlus.format(dt, "yyyy-e-mm-dd") == "2014-Sat-08-23"

    @test DatesPlus.format(DatesPlus.DateTime(2014, 1, 2, 0, 0, 0, 999), DatesPlus.RFC1123Format) == "Thu, 02 Jan 2014 00:00:00"
    @test DatesPlus.format(DatesPlus.DateTime(2014, 2, 18, 0, 0, 0, 9), DatesPlus.RFC1123Format) == "Tue, 18 Feb 2014 00:00:00"
    @test DatesPlus.format(DatesPlus.DateTime(2014, 3, 8, 0, 0, 0, 9), DatesPlus.RFC1123Format) == "Sat, 08 Mar 2014 00:00:00"
    @test DatesPlus.format(DatesPlus.DateTime(2014, 4, 28, 0, 0, 0, 9), DatesPlus.RFC1123Format) == "Mon, 28 Apr 2014 00:00:00"
    @test DatesPlus.format(DatesPlus.DateTime(2014, 5, 10, 0, 0, 0, 9), DatesPlus.RFC1123Format) == "Sat, 10 May 2014 00:00:00"
    @test DatesPlus.format(DatesPlus.DateTime(2014, 6, 4, 0, 0, 0, 9), DatesPlus.RFC1123Format) == "Wed, 04 Jun 2014 00:00:00"
    @test DatesPlus.format(DatesPlus.DateTime(2014, 7, 13, 0, 0, 0, 9), DatesPlus.RFC1123Format) == "Sun, 13 Jul 2014 00:00:00"
    @test DatesPlus.format(DatesPlus.DateTime(2014, 8, 17, 0, 0, 0, 9), DatesPlus.RFC1123Format) == "Sun, 17 Aug 2014 00:00:00"
    @test DatesPlus.format(DatesPlus.DateTime(2014, 9, 20, 0, 0, 0, 9), DatesPlus.RFC1123Format) == "Sat, 20 Sep 2014 00:00:00"
    @test DatesPlus.format(DatesPlus.DateTime(2014, 10, 31, 0, 0, 0, 9), DatesPlus.RFC1123Format) == "Fri, 31 Oct 2014 00:00:00"
    @test DatesPlus.format(DatesPlus.DateTime(2014, 11, 2, 0, 0, 0, 9), DatesPlus.RFC1123Format) == "Sun, 02 Nov 2014 00:00:00"
    @test DatesPlus.format(DatesPlus.DateTime(2014, 12, 5, 0, 0, 0, 9), DatesPlus.RFC1123Format) == "Fri, 05 Dec 2014 00:00:00"

    dt = DatesPlus.DateTime(2016, 11, 12, 7, 45, 36)
    @test parse(DatesPlus.DateTime, "Sat, 12 Nov 2016 07:45:36", DatesPlus.RFC1123Format) == dt
    @test parse(DatesPlus.DateTime, "Mon, 12 Nov 2016 07:45:36", DatesPlus.RFC1123Format) == dt  # Wrong day of week
    @test_throws ArgumentError parse(Date, "Foo, 12 Nov 2016 07:45:36", DatesPlus.RFC1123Format)
end

@testset "Issue 15195" begin
    f = "YY"
    @test DatesPlus.format(DatesPlus.Date(1999), f) == "1999"
    @test DatesPlus.format(DatesPlus.Date(9), f) == "09"
    @test DatesPlus.format(typemax(DatesPlus.Date), f) == "252522163911149"
end

@testset "TimeZones.jl Issue #19" begin
    Zulu = String

    function DatesPlus.tryparsenext(d::DatesPlus.DatePart{'Z'}, str, i, len)
        return DatesPlus.tryparsenext_word(str, i, len, DatesPlus.min_width(d), DatesPlus.max_width(d))
    end

    str = "2015-07-24T05:38:19.591Z"
    dt = DatesPlus.DateTime(2015, 7, 24, 5, 38, 19, 591)
    parsed = Any[
        DatesPlus.Year(2015), DatesPlus.Month(7), DatesPlus.Day(24),
        DatesPlus.Hour(5), DatesPlus.Minute(38), DatesPlus.Second(19), DatesPlus.Millisecond(591)
    ]

    format = "yyyy-mm-ddTHH:MM:SS.sssZ"
    escaped_format = "yyyy-mm-dd\\THH:MM:SS.sss\\Z"

    # Typically 'Z' isn't treated as a specifier so it doesn't have to be escaped
    @test DatesPlus.parse_components(str, DatesPlus.DateFormat(format)) == parsed
    @test DatesPlus.parse_components(str, DatesPlus.DateFormat(escaped_format)) == parsed

    try
        # Make 'Z' into a specifier
        DatesPlus.CONVERSION_SPECIFIERS['Z'] = Zulu
        DatesPlus.CONVERSION_DEFAULTS[Zulu] = ""

        @test DatesPlus.parse_components(str, DatesPlus.DateFormat(format)) == [parsed; Zulu("Z")]
        @test DatesPlus.parse_components(str, DatesPlus.DateFormat(escaped_format)) == parsed
    finally
        delete!(DatesPlus.CONVERSION_SPECIFIERS, 'Z')
        delete!(DatesPlus.CONVERSION_DEFAULTS, Zulu)
    end

    # Ensure that the default behaviour has been restored
    @test DatesPlus.parse_components(str, DatesPlus.DateFormat(format)) == parsed
    @test DatesPlus.parse_components(str, DatesPlus.DateFormat(escaped_format)) == parsed
end

@testset "Issue 10817" begin
    @test DatesPlus.Date("Apr 01 2014", "uuu dd yyyy") == DatesPlus.Date(2014, 4, 1)
    @test_throws ArgumentError DatesPlus.Date("Apr 01 xx 2014", "uuu dd zz yyyy")
    @test_throws ArgumentError DatesPlus.Date("Apr 01 xx 2014", "uuu dd    yyyy")
end

@testset "Issue 21001" begin
    for (ms, str) in zip([0, 1, 20, 300, 450, 678], ["0", "001", "02", "3", "45", "678"])
        local dt = DateTime(2000, 1, 1, 0, 0, 0, ms)
        @test DatesPlus.format(dt, "s") == str
        @test DatesPlus.format(dt, "ss") == rpad(str, 2, '0')
        @test DatesPlus.format(dt, "sss") == rpad(str, 3, '0')
        @test DatesPlus.format(dt, "ssss") == rpad(str, 4, '0')
    end
end

# Issue #21504
@test tryparse(DatesPlus.Date, "0-1000") === nothing

# Issue #44003
@test tryparse(DatesPlus.Date, "2017", DatesPlus.DateFormat(".s")) === nothing

@testset "parse milliseconds, Issue #22100" begin
    @test DatesPlus.DateTime("2017-Mar-17 00:00:00.0000", "y-u-d H:M:S.s") == DatesPlus.DateTime(2017, 3, 17)
    @test DatesPlus.parse_components(".1", DatesPlus.DateFormat(".s")) == [DatesPlus.Millisecond(100)]
    @test DatesPlus.parse_components(".12", DatesPlus.DateFormat(".s")) == [DatesPlus.Millisecond(120)]
    @test DatesPlus.parse_components(".123", DatesPlus.DateFormat(".s")) == [DatesPlus.Millisecond(123)]
    @test DatesPlus.parse_components(".1230", DatesPlus.DateFormat(".s")) == [DatesPlus.Millisecond(123)]
    # Issue #44003
    @test_throws ArgumentError DatesPlus.parse_components(".1234", DatesPlus.DateFormat(".s"))

    # Ensure that no overflow occurs when using Int32 literals: Int32(10)^10
    @test DatesPlus.parse_components("." * rpad(999, 10, '0'), DatesPlus.DateFormat(".s")) == [DatesPlus.Millisecond(999)]
end

@testset "Time Parsing" begin
    let t
        time_tuple(t::DatesPlus.Time) = (
            DatesPlus.hour(t), DatesPlus.minute(t), DatesPlus.second(t),
            DatesPlus.millisecond(t), DatesPlus.microsecond(t), DatesPlus.nanosecond(t)
        )

        ## default ISOTimeFormat
        t = DatesPlus.Time("01")
        @test time_tuple(t) == (1, 0, 0, 0, 0, 0)
        t = DatesPlus.Time("01:23")
        @test time_tuple(t) == (1, 23, 0, 0, 0, 0)
        t = DatesPlus.Time("01:23:45")
        @test time_tuple(t) == (1, 23, 45, 0, 0, 0)
        t = DatesPlus.Time("01:23:45.678")
        @test time_tuple(t) == (1, 23, 45, 678, 0, 0)

        ## string format
        t = DatesPlus.Time("23:56:12.1", "HH:MM:SS.s")
        @test time_tuple(t) == (23, 56, 12, 100, 0, 0)

        ## precomputed DateFormat
        t = DatesPlus.Time("04:09:45.012", DateFormat("HH:MM:SS.s"))
        @test time_tuple(t) == (4, 9, 45, 12, 0, 0)
        t = DatesPlus.Time("21 07", DateFormat("HH MM"))
        @test time_tuple(t) == (21, 7, 0, 0, 0, 0)
        t = DatesPlus.Time("4.02", DateFormat("H.MM"))
        @test time_tuple(t) == (4, 2, 0, 0, 0, 0)
        t = DatesPlus.Time("1725", DateFormat("HHMM"))
        @test time_tuple(t) == (17, 25, 0, 0, 0, 0)

        ## exceptions
        @test_throws ArgumentError DatesPlus.Time("24:00")  # invalid hours
        @test_throws ArgumentError DatesPlus.Time("00:60")  # invalid minutes
        @test_throws ArgumentError DatesPlus.Time("00:00:60")  # invalid seconds
        @test_throws ArgumentError DatesPlus.Time("20:03:20", DateFormat("HH:MM"))  # too much precision
        @test_throws ArgumentError DatesPlus.Time("10:33:51", DateFormat("YYYY-MM-DD HH:MM:SS"))  # Time can't hold year/month/day
    end
end

@testset "midnight" begin
    # issue #28203: 24:00 is a valid ISO 8601 time
    @test DateTime("2018-01-01 24:00", "yyyy-mm-dd HH:MM") ==
          DateTime("2018-01-02T00:00:00") ==
          DateTime(2018, 1, 1, 24) ==
          DateTime(2018, 1, 2)
    @test_throws ArgumentError DateTime("2018-01-01 24:01", "yyyy-mm-dd HH:MM")
    @test_throws ArgumentError DateTime(2018, 1, 1, 24, 0, 1)
    @test_throws ArgumentError DateTime(2018, 1, 1, 24, 0, 0, 1)
end

@testset "AM/PM" begin
    for (t12, t24) in (
        ("12:00am", "00:00"),
        ("12:07am", "00:07"),
        ("01:24AM", "01:24"),
        ("12:00pm", "12:00"),
        ("12:15pm", "12:15"),
        ("11:59PM", "23:59"),
    )
        d = DateTime("2018-01-01T$t24:00")
        t = Time("$t24:00")
        for HH in ("HH", "II")
            @test DateTime("2018-01-01 $t12", "yyyy-mm-dd $HH:MMp") == d
            @test Time("$t12", "$HH:MMp") == t
        end
        local tmstruct, strftime
        withlocales(["C"]) do
            # test am/pm comparison handling
            tmstruct = Libc.strptime("%I:%M%p", t12)
            strftime = Libc.strftime("%I:%M%p", tmstruct)
            nothing
        end
        @test Time(tmstruct) == t
        @test uppercase(t12) ==
              DatesPlus.format(t, "II:MMp") ==
              DatesPlus.format(d, "II:MMp") ==
              strftime
    end
    for bad in ("00:24am", "00:24pm", "13:24pm", "2pm", "12:24p.m.", "12:24 pm", "12:24pµ")
        @test_throws ArgumentError Time(bad, "II:MMp")
    end
    # if am/pm is missing, defaults to 24-hour clock
    @test Time("13:24", "II:MMp") == Time("13:24", "HH:MM")
end

@testset "Issue #10561, two-digit year parsing ambiguities" begin
    # All two-digit dates (whether full or resulting from truncation, e.g. 2010 -> 10)
    # encoded in two digit year format YY are parsed as year 00YY

    for test_year in ("00", "01", "99", "2021")
        @test Date(test_year, dateformat"yy") == Date(parse(Int, test_year))
    end
end

@testset "inference with dynamic dateformat string" begin
    datetime = DateTime(2020, 4, 7)
    f1() = DateTime("2020-04-07", "yyyy-mm-dd")
    f2() = DateTime("2020-04-07", DateFormat("yyyy-mm-dd"))
    f3() = parse(DateTime, "2020-04-07", DateFormat("yyyy-mm-dd"))
    @test (@inferred f1()) == (@inferred f2()) == (@inferred f3()) == datetime
    g() = tryparse(DateTime, "2020-04-07", DateFormat("yyyy-mm-dd"))
    @test (@inferred Nothing g()) == datetime
end

end


text:
"""
    Collection(ref::Vector{Reference})

Wrapper type around all kind of Reference types. Is used to define special methods for vectors of references.
"""
mutable struct Collection <: Reference
    refs::Vector{Reference}
end

==(a::Collection,b::Collection) = a.refs == b.refs

@alive c value(c::Collection) = value.(c.refs)

@alive c function destroy!(c::Collection)
    destroy!.(c.refs)
    c.refs = Vector{Reference}()
end
isdestroyed(c::Collection) = length(c.refs) == 0

function show(io::IO,c::Collection)
    if isdestroyed(c)
        print(io,"*Deleted*")
    else
        print(io,"$(length(c.refs))-element Collection:")
        for l in c.refs
            print(io,"\n $l")
        end
    end
end

iterate(c::Collection) = (c[1],1)

function iterate(c::Collection,state::Int)
    if state < length(c)
        return (c[state+1],state+1)
    else
        return nothing
    end
end

getindex(c::Collection,i::Int) = c.refs[i]

length(c::Collection) = length(c.refs)

@alive c next(c::Collection,i::Int=1) = next.(c.refs,i)
@alive c prev(c::Collection,i::Int=1) = prev.(c.refs,i)

+(c::Collection,i::Int) = next(c,i)
-(c::Collection,i::Int) = prev(c,i)
+(i::Int,c::Collection) = next(c,i)
-(i::Int,c::Collection) = prev(c,i)

@alive c copy(c::Collection) = Collection(copy(c.refs))

@alive c delete!(c::Collection) = delete!.(c.refs)
notify_delete(c::Collection,del::T) where T<:Reference = notify_delete.(c.refs,Ref{T}(del))

@alive c from insert!(c::Collection,from::T) where T<:Reference = insert!.(c.refs,Ref{T}(from)) |> Collection
@alive c insert!(c::Collection,from::T) where T<:Union{String,Vector{String}} = insert!.(c.refs,Ref{T}(from)) |> Collection
notify_insert(self::Collection,ins::T,lines::Int) where T<:Reference = notify_insert.(self.refs,Ref{T}(ins),Ref{T}(lines))

@alive c from append!(c::Collection,from::T) where T<:Reference = append!.(c.refs,from) |> Collection
@alive c append!(c::Collection,from::T) where T<:Union{String,Vector{String}} = append!.(c.refs,Ref{T}(from)) |> Collection
notify_append(self::Collection,app::T,lines::Int) where T<:Reference = notify_append.(self.refs,Ref{T}(app),lines)

@alive c from replace!(c::Collection,from::T) where T<:Reference = replace!.(c.refs,from) |> Collection
@alive c replace!(c::Collection,from::T) where T<:Union{String,Vector{String}} = replace!.(c.refs,Ref{T}(from)) |> Collection
notify_replace(self::Collection,rep::T,lines::Int) where T<:Reference = notify_replace.(self.refs,Ref{T}(rep),lines)

@alive self to function move!(self::Collection,to::Collection)
    length(self.refs) != length(to.refs) && error("Collection can only be moved to other collections with equal length.")
    for i in 1:length(self.refs)
        move!(self.refs[i],to.refs[i])
    end
end

@alive self to function moveafter!(self::Collection,to::Collection)
    length(self.refs) != length(to.refs) && error("Collection can only be moved to other collections with equal length.")
    for i in 1:length(self.refs)
        moveafter!(self.refs[i],to.refs[i])
    end
end


text:
export  LinHermiteMapComponent,
        getcoeff,
        setcoeff!,
        getidx,
        evaluate!,
        evaluate

"""
$(TYPEDEF)

`LinHermiteMapComponent` is a type to hold  the composition of a `LinearTransform`: a linear transformation that rescales
the samples to get zero mean and unitary standard deviation component-wise
and  an `HermiteMapComponent` element.

## Fields

$(TYPEDFIELDS)
"""
struct LinHermiteMapComponent
    # Linear transformation
    L::LinearTransform

    # IntegratedFunction
    C::HermiteMapComponent

    function LinHermiteMapComponent(L::LinearTransform, C::HermiteMapComponent)
        return new(L, C)
    end
end

function LinHermiteMapComponent(X::Array{Float64,2}, C::HermiteMapComponent; diag::Bool=true)
    @assert size(X,1)==C.Nx "Wrong dimension of the input"
    L = LinearTransform(X; diag = diag)
    return LinHermiteMapComponent(L, C)
end

"""
$(TYPEDSIGNATURES)

Returns the number of multi-variate features of `L`.
"""
ncoeff(L::LinHermiteMapComponent) = L.C.Nψ
"""
$(TYPEDSIGNATURES)

Returns the coefficients of `L`.
"""
getcoeff(L::LinHermiteMapComponent) = L.C.I.f.coeff

"""
$(TYPEDSIGNATURES)

Sets the coefficients of `L`  to `coeff`.
"""
function setcoeff!(L::LinHermiteMapComponent, coeff::Array{Float64,1})
        @assert size(coeff,1) == Nψ "Wrong dimension of coeff"
        L.C.I.f.coeff .= coeff
end

"""
$(TYPEDSIGNATURES)

Extracts the multi-indices of `L`.
"""
getidx(L::LinHermiteMapComponent) = L.C.I.f.idx

"""
$(TYPEDSIGNATURES)

Evaluates in-place the `LinHermiteMapComponent` `L` for the ensemble matrix `X`.
"""
function evaluate!(out, L::LinHermiteMapComponent, X)
    @assert L.C.Nx==size(X,1) "Wrong dimension of the sample"
    @assert size(out,1) == size(X,2) "Dimensions of the output and the samples don't match"

    transform!(L.L, X)
    evaluate!(out, L.C.I, X)
    itransform!(L.L, X)
    return out
end

"""
$(TYPEDSIGNATURES)

Evaluates the `LinHermiteMapComponent` `L` for the ensemble matrix `X`.
"""
evaluate(L::LinHermiteMapComponent, X::Array{Float64,2}) = evaluate!(zeros(size(X,2)), L, X)


text:
function f(v :: String, w :: Int = 12, x :: Number, y = 11, z ; a = "chae", b = 2)
  x + y
end

max_stars_repo_name: dadarakt/jamon

max_issues_repo_name: dadarakt/jamon

max_forks_repo_name: dadarakt/jamon


text:
# This problem gives practice on reusing procedures by extracting them as functions.

function is_equilateral(sides)
    sides = sort(sides)
    return verify_triangle(sides) && verify_equilateral(sides)
end

function is_isosceles(sides)
    sides = sort(sides)
    return verify_triangle(sides) && verify_isosceles(sides)
end

function is_scalene(sides)
    sides = sort(sides)
    return verify_triangle(sides) && verify_scalene(sides)
end

verify_equilateral(sorted) = (sorted[1] == sorted[3])

verify_isosceles(sorted) = (sorted[1] == sorted[2] || sorted[2] == sorted[3])

verify_scalene(sorted) = !verify_isosceles(sorted)

verify_triangle(sorted) = (sorted[3] < sorted[1] + sorted[2])


text:
# graphstructure.jl
# Authors : Heitor Guimaraes and Luiz Ciafrino
# @brief: Module to handle basic plot for the graph algorithm's

using Colors

function plot_graph(G::simpleGraph)
	println(" -- Shape of the Graph --")
end

function plot_bfs()
	println("-- Result of the BFS --")
end

function plot_dfs()
	println(" -- Result of the DFS --")
end


text:
#!/bin/env julia
import Pkg
Pkg.activate("./Threaded")

Enzyme_path = Base.find_package("Enzyme")
if Enzyme_path === nothing
    error("First instantiate Project")
end
dir = joinpath(dirname(dirname(Enzyme_path)), "deps")

@info "Building Enzyme from inc"

run(`$(Base.julia_cmd()) --project=$(dir) -e 'import Pkg; Pkg.instantiate()'`)
run(`$(Base.julia_cmd()) --project=$(dir) $(dir)/build_local.jl --branch vc/jl_inc`)

cp(joinpath(dirname(dir), "LocalPreferences.toml"), "Threaded/LocalPreferences.toml", force=true)

max_stars_repo_name: wsmoses/Enzyme-BUDE

max_issues_repo_name: wsmoses/Enzyme-BUDE

max_forks_repo_name: wsmoses/Enzyme-BUDE