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

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

Creates a broker. Note: This API is asynchronous.

# Optional Parameters
Optional parameters can be passed as a `params::Dict{String,<:Any}`. Valid keys are:
- `"authenticationStrategy"`: The authentication strategy used to secure the broker.
- `"autoMinorVersionUpgrade"`: Required. Enables automatic upgrades to new minor versions
  for brokers, as Apache releases the versions. The automatic upgrades occur during the
  maintenance window of the broker or after a manual broker reboot.
- `"brokerName"`: Required. The name of the broker. This value must be unique in your AWS
  account, 1-50 characters long, must contain only letters, numbers, dashes, and underscores,
  and must not contain whitespaces, brackets, wildcard characters, or special characters.
- `"configuration"`: A list of information about the configuration.
- `"creatorRequestId"`: The unique ID that the requester receives for the created broker.
  Amazon MQ passes your ID with the API action. Note: We recommend using a Universally Unique
  Identifier (UUID) for the creatorRequestId. You may omit the creatorRequestId if your
  application doesn't require idempotency.
- `"deploymentMode"`: Required. The deployment mode of the broker.
- `"encryptionOptions"`: Encryption options for the broker.
- `"engineType"`: Required. The type of broker engine. Note: Currently, Amazon MQ supports
  ACTIVEMQ and RABBITMQ.
- `"engineVersion"`: Required. The version of the broker engine. For a list of supported
  engine versions, see
  https://docs.aws.amazon.com/amazon-mq/latest/developer-guide/broker-engine.html
- `"hostInstanceType"`: Required. The broker's instance type.
- `"ldapServerMetadata"`: The metadata of the LDAP server used to authenticate and
  authorize connections to the broker.
- `"logs"`: Enables Amazon CloudWatch logging for brokers.
- `"maintenanceWindowStartTime"`: The parameters that determine the WeeklyStartTime.
- `"publiclyAccessible"`: Required. Enables connections from applications outside of the
  VPC that hosts the broker's subnets.
- `"securityGroups"`: The list of security groups (1 minimum, 5 maximum) that authorizes
  connections to brokers.
- `"storageType"`: The broker's storage type.
- `"subnetIds"`: The list of groups that define which subnets and IP ranges the broker can
  use from different Availability Zones. A SINGLE_INSTANCE deployment requires one subnet
  (for example, the default subnet). An ACTIVE_STANDBY_MULTI_AZ deployment (ACTIVEMQ)
  requires two subnets. A CLUSTER_MULTI_AZ deployment (RABBITMQ) has no subnet requirements
  when deployed with public accessibility, deployment without public accessibility requires
  at least one subnet.
- `"tags"`: Create tags when creating the broker.
- `"users"`: Required. The list of broker users (persons or applications) who can access
  queues and topics. For RabbitMQ brokers, one and only one administrative user is accepted
  and created when a broker is first provisioned. All subsequent broker users are created by
  making RabbitMQ API calls directly to brokers or via the RabbitMQ Web Console. This value
  can contain only alphanumeric characters, dashes, periods, underscores, and tildes (- . _
  ~). This value must be 2-100 characters long.
"""
create_broker(; aws_config::AbstractAWSConfig=global_aws_config()) = mq("POST", "/v1/brokers", Dict{String, Any}("creatorRequestId"=>string(uuid4())); aws_config=aws_config)
create_broker(params::AbstractDict{String}; aws_config::AbstractAWSConfig=global_aws_config()) = mq("POST", "/v1/brokers", Dict{String, Any}(mergewith(_merge, Dict{String, Any}("creatorRequestId"=>string(uuid4())), params)); aws_config=aws_config)

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

Creates a new configuration for the specified configuration name. Amazon MQ uses the
default configuration (the engine type and version).

# Optional Parameters
Optional parameters can be passed as a `params::Dict{String,<:Any}`. Valid keys are:
- `"authenticationStrategy"`: The authentication strategy associated with the configuration.
- `"engineType"`: Required. The type of broker engine. Note: Currently, Amazon MQ supports
  ACTIVEMQ and RABBITMQ.
- `"engineVersion"`: Required. The version of the broker engine. For a list of supported
  engine versions, see
  https://docs.aws.amazon.com/amazon-mq/latest/developer-guide/broker-engine.html
- `"name"`: Required. The name of the configuration. This value can contain only
  alphanumeric characters, dashes, periods, underscores, and tildes (- . _ ~). This value
  must be 1-150 characters long.
- `"tags"`: Create tags when creating the configuration.
"""
create_configuration(; aws_config::AbstractAWSConfig=global_aws_config()) = mq("POST", "/v1/configurations"; aws_config=aws_config)
create_configuration(params::AbstractDict{String}; aws_config::AbstractAWSConfig=global_aws_config()) = mq("POST", "/v1/configurations", params; aws_config=aws_config)

"""
    create_tags(resource-arn)
    create_tags(resource-arn, params::Dict{String,<:Any})

Add a tag to a resource.

# Arguments
- `resource-arn`: The Amazon Resource Name (ARN) of the resource tag.

# Optional Parameters
Optional parameters can be passed as a `params::Dict{String,<:Any}`. Valid keys are:
- `"tags"`: The key-value pair for the resource tag.
"""
create_tags(resource_arn; aws_config::AbstractAWSConfig=global_aws_config()) = mq("POST", "/v1/tags/$(resource-arn)"; aws_config=aws_config)
create_tags(resource_arn, params::AbstractDict{String}; aws_config::AbstractAWSConfig=global_aws_config()) = mq("POST", "/v1/tags/$(resource-arn)", params; aws_config=aws_config)

"""
    create_user(broker-id, username)
    create_user(broker-id, username, params::Dict{String,<:Any})

Creates an ActiveMQ user.

# Arguments
- `broker-id`: The unique ID that Amazon MQ generates for the broker.
- `username`: The username of the ActiveMQ user. This value can contain only alphanumeric
  characters, dashes, periods, underscores, and tildes (- . _ ~). This value must be 2-100
  characters long.

# Optional Parameters
Optional parameters can be passed as a `params::Dict{String,<:Any}`. Valid keys are:
- `"consoleAccess"`: Enables access to the the ActiveMQ Web Console for the ActiveMQ user.
- `"groups"`: The list of groups (20 maximum) to which the ActiveMQ user belongs. This
  value can contain only alphanumeric characters, dashes, periods, underscores, and tildes (-
  . _ ~). This value must be 2-100 characters long.
- `"password"`: Required. The password of the user. This value must be at least 12
  characters long, must contain at least 4 unique characters, and must not contain commas.
"""
create_user(broker_id, username; aws_config::AbstractAWSConfig=global_aws_config()) = mq("POST", "/v1/brokers/$(broker-id)/users/$(username)"; aws_config=aws_config)
create_user(broker_id, username, params::AbstractDict{String}; aws_config::AbstractAWSConfig=global_aws_config()) = mq("POST", "/v1/brokers/$(broker-id)/users/$(username)", params; aws_config=aws_config)

"""
    delete_broker(broker-id)
    delete_broker(broker-id, params::Dict{String,<:Any})

Deletes a broker. Note: This API is asynchronous.

# Arguments
- `broker-id`: The unique ID that Amazon MQ generates for the broker.

"""
delete_broker(broker_id; aws_config::AbstractAWSConfig=global_aws_config()) = mq("DELETE", "/v1/brokers/$(broker-id)"; aws_config=aws_config)
delete_broker(broker_id, params::AbstractDict{String}; aws_config::AbstractAWSConfig=global_aws_config()) = mq("DELETE", "/v1/brokers/$(broker-id)", params; aws_config=aws_config)

"""
    delete_tags(resource-arn, tag_keys)
    delete_tags(resource-arn, tag_keys, params::Dict{String,<:Any})

Removes a tag from a resource.

# Arguments
- `resource-arn`: The Amazon Resource Name (ARN) of the resource tag.
- `tag_keys`: An array of tag keys to delete

"""
delete_tags(resource_arn, tagKeys; aws_config::AbstractAWSConfig=global_aws_config()) = mq("DELETE", "/v1/tags/$(resource-arn)", Dict{String, Any}("tagKeys"=>tagKeys); aws_config=aws_config)
delete_tags(resource_arn, tagKeys, params::AbstractDict{String}; aws_config::AbstractAWSConfig=global_aws_config()) = mq("DELETE", "/v1/tags/$(resource-arn)", Dict{String, Any}(mergewith(_merge, Dict{String, Any}("tagKeys"=>tagKeys), params)); aws_config=aws_config)

"""
    delete_user(broker-id, username)
    delete_user(broker-id, username, params::Dict{String,<:Any})

Deletes an ActiveMQ user.

# Arguments
- `broker-id`: The unique ID that Amazon MQ generates for the broker.
- `username`: The username of the ActiveMQ user. This value can contain only alphanumeric
  characters, dashes, periods, underscores, and tildes (- . _ ~). This value must be 2-100
  characters long.

"""
delete_user(broker_id, username; aws_config::AbstractAWSConfig=global_aws_config()) = mq("DELETE", "/v1/brokers/$(broker-id)/users/$(username)"; aws_config=aws_config)
delete_user(broker_id, username, params::AbstractDict{String}; aws_config::AbstractAWSConfig=global_aws_config()) = mq("DELETE", "/v1/brokers/$(broker-id)/users/$(username)", params; aws_config=aws_config)

"""
    describe_broker(broker-id)
    describe_broker(broker-id, params::Dict{String,<:Any})

Returns information about the specified broker.

# Arguments
- `broker-id`: The name of the broker. This value must be unique in your AWS account, 1-50
  characters long, must contain only letters, numbers, dashes, and underscores, and must not
  contain whitespaces, brackets, wildcard characters, or special characters.

"""
describe_broker(broker_id; aws_config::AbstractAWSConfig=global_aws_config()) = mq("GET", "/v1/brokers/$(broker-id)"; aws_config=aws_config)
describe_broker(broker_id, params::AbstractDict{String}; aws_config::AbstractAWSConfig=global_aws_config()) = mq("GET", "/v1/brokers/$(broker-id)", params; aws_config=aws_config)

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

Describe available engine types and versions.

# Optional Parameters
Optional parameters can be passed as a `params::Dict{String,<:Any}`. Valid keys are:
- `"engineType"`: Filter response by engine type.
- `"maxResults"`: The maximum number of engine types that Amazon MQ can return per page (20
  by default). This value must be an integer from 5 to 100.
- `"nextToken"`: The token that specifies the next page of results Amazon MQ should return.
  To request the first page, leave nextToken empty.
"""
describe_broker_engine_types(; aws_config::AbstractAWSConfig=global_aws_config()) = mq("GET", "/v1/broker-engine-types"; aws_config=aws_config)
describe_broker_engine_types(params::AbstractDict{String}; aws_config::AbstractAWSConfig=global_aws_config()) = mq("GET", "/v1/broker-engine-types", params; aws_config=aws_config)

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

Describe available broker instance options.

# Optional Parameters
Optional parameters can be passed as a `params::Dict{String,<:Any}`. Valid keys are:
- `"engineType"`: Filter response by engine type.
- `"hostInstanceType"`: Filter response by host instance type.
- `"maxResults"`: The maximum number of instance options that Amazon MQ can return per page
  (20 by default). This value must be an integer from 5 to 100.
- `"nextToken"`: The token that specifies the next page of results Amazon MQ should return.
  To request the first page, leave nextToken empty.
- `"storageType"`: Filter response by storage type.
"""
describe_broker_instance_options(; aws_config::AbstractAWSConfig=global_aws_config()) = mq("GET", "/v1/broker-instance-options"; aws_config=aws_config)
describe_broker_instance_options(params::AbstractDict{String}; aws_config::AbstractAWSConfig=global_aws_config()) = mq("GET", "/v1/broker-instance-options", params; aws_config=aws_config)

"""
    describe_configuration(configuration-id)
    describe_configuration(configuration-id, params::Dict{String,<:Any})

Returns information about the specified configuration.

# Arguments
- `configuration-id`: The unique ID that Amazon MQ generates for the configuration.

"""
describe_configuration(configuration_id; aws_config::AbstractAWSConfig=global_aws_config()) = mq("GET", "/v1/configurations/$(configuration-id)"; aws_config=aws_config)
describe_configuration(configuration_id, params::AbstractDict{String}; aws_config::AbstractAWSConfig=global_aws_config()) = mq("GET", "/v1/configurations/$(configuration-id)", params; aws_config=aws_config)

"""
    describe_configuration_revision(configuration-id, configuration-revision)
    describe_configuration_revision(configuration-id, configuration-revision, params::Dict{String,<:Any})

Returns the specified configuration revision for the specified configuration.

# Arguments
- `configuration-id`: The unique ID that Amazon MQ generates for the configuration.
- `configuration-revision`: The revision of the configuration.

"""
describe_configuration_revision(configuration_id, configuration_revision; aws_config::AbstractAWSConfig=global_aws_config()) = mq("GET", "/v1/configurations/$(configuration-id)/revisions/$(configuration-revision)"; aws_config=aws_config)
describe_configuration_revision(configuration_id, configuration_revision, params::AbstractDict{String}; aws_config::AbstractAWSConfig=global_aws_config()) = mq("GET", "/v1/configurations/$(configuration-id)/revisions/$(configuration-revision)", params; aws_config=aws_config)

"""
    describe_user(broker-id, username)
    describe_user(broker-id, username, params::Dict{String,<:Any})

Returns information about an ActiveMQ user.

# Arguments
- `broker-id`: The unique ID that Amazon MQ generates for the broker.
- `username`: The username of the ActiveMQ user. This value can contain only alphanumeric
  characters, dashes, periods, underscores, and tildes (- . _ ~). This value must be 2-100
  characters long.

"""
describe_user(broker_id, username; aws_config::AbstractAWSConfig=global_aws_config()) = mq("GET", "/v1/brokers/$(broker-id)/users/$(username)"; aws_config=aws_config)
describe_user(broker_id, username, params::AbstractDict{String}; aws_config::AbstractAWSConfig=global_aws_config()) = mq("GET", "/v1/brokers/$(broker-id)/users/$(username)", params; aws_config=aws_config)

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

Returns a list of all brokers.

# Optional Parameters
Optional parameters can be passed as a `params::Dict{String,<:Any}`. Valid keys are:
- `"maxResults"`: The maximum number of brokers that Amazon MQ can return per page (20 by
  default). This value must be an integer from 5 to 100.
- `"nextToken"`: The token that specifies the next page of results Amazon MQ should return.
  To request the first page, leave nextToken empty.
"""
list_brokers(; aws_config::AbstractAWSConfig=global_aws_config()) = mq("GET", "/v1/brokers"; aws_config=aws_config)
list_brokers(params::AbstractDict{String}; aws_config::AbstractAWSConfig=global_aws_config()) = mq("GET", "/v1/brokers", params; aws_config=aws_config)

"""
    list_configuration_revisions(configuration-id)
    list_configuration_revisions(configuration-id, params::Dict{String,<:Any})

Returns a list of all revisions for the specified configuration.

# Arguments
- `configuration-id`: The unique ID that Amazon MQ generates for the configuration.

# Optional Parameters
Optional parameters can be passed as a `params::Dict{String,<:Any}`. Valid keys are:
- `"maxResults"`: The maximum number of configurations that Amazon MQ can return per page
  (20 by default). This value must be an integer from 5 to 100.
- `"nextToken"`: The token that specifies the next page of results Amazon MQ should return.
  To request the first page, leave nextToken empty.
"""
list_configuration_revisions(configuration_id; aws_config::AbstractAWSConfig=global_aws_config()) = mq("GET", "/v1/configurations/$(configuration-id)/revisions"; aws_config=aws_config)
list_configuration_revisions(configuration_id, params::AbstractDict{String}; aws_config::AbstractAWSConfig=global_aws_config()) = mq("GET", "/v1/configurations/$(configuration-id)/revisions", params; aws_config=aws_config)

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

Returns a list of all configurations.

# Optional Parameters
Optional parameters can be passed as a `params::Dict{String,<:Any}`. Valid keys are:
- `"maxResults"`: The maximum number of configurations that Amazon MQ can return per page
  (20 by default). This value must be an integer from 5 to 100.
- `"nextToken"`: The token that specifies the next page of results Amazon MQ should return.
  To request the first page, leave nextToken empty.
"""
list_configurations(; aws_config::AbstractAWSConfig=global_aws_config()) = mq("GET", "/v1/configurations"; aws_config=aws_config)
list_configurations(params::AbstractDict{String}; aws_config::AbstractAWSConfig=global_aws_config()) = mq("GET", "/v1/configurations", params; aws_config=aws_config)

"""
    list_tags(resource-arn)
    list_tags(resource-arn, params::Dict{String,<:Any})

Lists tags for a resource.

# Arguments
- `resource-arn`: The Amazon Resource Name (ARN) of the resource tag.

"""
list_tags(resource_arn; aws_config::AbstractAWSConfig=global_aws_config()) = mq("GET", "/v1/tags/$(resource-arn)"; aws_config=aws_config)
list_tags(resource_arn, params::AbstractDict{String}; aws_config::AbstractAWSConfig=global_aws_config()) = mq("GET", "/v1/tags/$(resource-arn)", params; aws_config=aws_config)

"""
    list_users(broker-id)
    list_users(broker-id, params::Dict{String,<:Any})

Returns a list of all ActiveMQ users.

# Arguments
- `broker-id`: The unique ID that Amazon MQ generates for the broker.

# Optional Parameters
Optional parameters can be passed as a `params::Dict{String,<:Any}`. Valid keys are:
- `"maxResults"`: The maximum number of ActiveMQ users that can be returned per page (20 by
  default). This value must be an integer from 5 to 100.
- `"nextToken"`: The token that specifies the next page of results Amazon MQ should return.
  To request the first page, leave nextToken empty.
"""
list_users(broker_id; aws_config::AbstractAWSConfig=global_aws_config()) = mq("GET", "/v1/brokers/$(broker-id)/users"; aws_config=aws_config)
list_users(broker_id, params::AbstractDict{String}; aws_config::AbstractAWSConfig=global_aws_config()) = mq("GET", "/v1/brokers/$(broker-id)/users", params; aws_config=aws_config)

"""
    reboot_broker(broker-id)
    reboot_broker(broker-id, params::Dict{String,<:Any})

Reboots a broker. Note: This API is asynchronous.

# Arguments
- `broker-id`: The unique ID that Amazon MQ generates for the broker.

"""
reboot_broker(broker_id; aws_config::AbstractAWSConfig=global_aws_config()) = mq("POST", "/v1/brokers/$(broker-id)/reboot"; aws_config=aws_config)
reboot_broker(broker_id, params::AbstractDict{String}; aws_config::AbstractAWSConfig=global_aws_config()) = mq("POST", "/v1/brokers/$(broker-id)/reboot", params; aws_config=aws_config)

"""
    update_broker(broker-id)
    update_broker(broker-id, params::Dict{String,<:Any})

Adds a pending configuration change to a broker.

# Arguments
- `broker-id`: The unique ID that Amazon MQ generates for the broker.

# Optional Parameters
Optional parameters can be passed as a `params::Dict{String,<:Any}`. Valid keys are:
- `"authenticationStrategy"`: The authentication strategy used to secure the broker.
- `"autoMinorVersionUpgrade"`: Enables automatic upgrades to new minor versions for
  brokers, as Apache releases the versions. The automatic upgrades occur during the
  maintenance window of the broker or after a manual broker reboot.
- `"configuration"`: A list of information about the configuration.
- `"engineVersion"`: The version of the broker engine. For a list of supported engine
  versions, see
  https://docs.aws.amazon.com/amazon-mq/latest/developer-guide/broker-engine.html
- `"hostInstanceType"`: The host instance type of the broker to upgrade to. For a list of
  supported instance types, see
  https://docs.aws.amazon.com/amazon-mq/latest/developer-guide//broker.html#broker-instance-ty
  pes
- `"ldapServerMetadata"`: The metadata of the LDAP server used to authenticate and
  authorize connections to the broker.
- `"logs"`: Enables Amazon CloudWatch logging for brokers.
- `"securityGroups"`: The list of security groups (1 minimum, 5 maximum) that authorizes
  connections to brokers.
"""
update_broker(broker_id; aws_config::AbstractAWSConfig=global_aws_config()) = mq("PUT", "/v1/brokers/$(broker-id)"; aws_config=aws_config)
update_broker(broker_id, params::AbstractDict{String}; aws_config::AbstractAWSConfig=global_aws_config()) = mq("PUT", "/v1/brokers/$(broker-id)", params; aws_config=aws_config)

"""
    update_configuration(configuration-id)
    update_configuration(configuration-id, params::Dict{String,<:Any})

Updates the specified configuration.

# Arguments
- `configuration-id`: The unique ID that Amazon MQ generates for the configuration.

# Optional Parameters
Optional parameters can be passed as a `params::Dict{String,<:Any}`. Valid keys are:
- `"data"`: Required. The base64-encoded XML configuration.
- `"description"`: The description of the configuration.
"""
update_configuration(configuration_id; aws_config::AbstractAWSConfig=global_aws_config()) = mq("PUT", "/v1/configurations/$(configuration-id)"; aws_config=aws_config)
update_configuration(configuration_id, params::AbstractDict{String}; aws_config::AbstractAWSConfig=global_aws_config()) = mq("PUT", "/v1/configurations/$(configuration-id)", params; aws_config=aws_config)

"""
    update_user(broker-id, username)
    update_user(broker-id, username, params::Dict{String,<:Any})

Updates the information for an ActiveMQ user.

# Arguments
- `broker-id`: The unique ID that Amazon MQ generates for the broker.
- `username`: Required. The username of the ActiveMQ user. This value can contain only
  alphanumeric characters, dashes, periods, underscores, and tildes (- . _ ~). This value
  must be 2-100 characters long.

# Optional Parameters
Optional parameters can be passed as a `params::Dict{String,<:Any}`. Valid keys are:
- `"consoleAccess"`: Enables access to the the ActiveMQ Web Console for the ActiveMQ user.
- `"groups"`: The list of groups (20 maximum) to which the ActiveMQ user belongs. This
  value can contain only alphanumeric characters, dashes, periods, underscores, and tildes (-
  . _ ~). This value must be 2-100 characters long.
- `"password"`: The password of the user. This value must be at least 12 characters long,
  must contain at least 4 unique characters, and must not contain commas.
"""
update_user(broker_id, username; aws_config::AbstractAWSConfig=global_aws_config()) = mq("PUT", "/v1/brokers/$(broker-id)/users/$(username)"; aws_config=aws_config)
update_user(broker_id, username, params::AbstractDict{String}; aws_config::AbstractAWSConfig=global_aws_config()) = mq("PUT", "/v1/brokers/$(broker-id)/users/$(username)", params; aws_config=aws_config)


text:
# `BENCHMARKS`, `suite`, and `add_benchmark` were taken from MathOptInterface
# which is available under an MIT license (see LICENSE).
module Benchmarks

using BenchmarkTools
using Convex, LinearAlgebra

const BENCHMARKS = Dict{String, Function}()

"""
    suite(
        handle_problem::Function;
        exclude::Vector{Regex} = Regex[]
    )

Create a suite of benchmarks. `handle_problem` should be a function that takes one
argument, a Convex.jl `Problem` and processes it (e.g. `solve!` the problem with
a specific solver).

Use `exclude` to exclude a subset of benchmarks.

### Examples

```julia
suite() do p
    solve!(p, GLPK.Optimizer())
end
```
"""
function suite(handle_problem::Function; exclude::Vector{Regex} = Regex[])
    group = BenchmarkGroup()
    for (name, func) in BENCHMARKS
        any(occursin.(exclude, Ref(name))) && continue
        group[name] = @benchmarkable $func($handle_problem) setup=Convex.clearmemory()
    end
    return group
end

###
### Benchmarks
###

macro add_benchmark(f)
    name = f.args[1].args[1]
    return quote
        $(esc(f))
        BENCHMARKS[String($(Base.Meta.quot(name)))] = $(esc(name))
    end
end

eye(n) = Matrix(1.0I, n, n)
include("constraint_benchmarks.jl")
include("affine_benchmarks.jl")
include("sdp_benchmarks.jl")

end

ext: jl


text:
module WaveletUtils

using StaticArrays

include("collections.jl")
include("detail.jl")
include("trees.jl")

end # module


text:
__precompile__(true)
module WeakRefStrings

export WeakRefString, WeakRefStringArray

using Missings

if !isdefined(Base, :codeunits)
    codeunits = Vector{UInt8}
end

"""
A custom "weakref" string type that only points to external string data.
Allows for the creation of a "string" instance without copying data,
which allows for more efficient string parsing/movement in certain data processing tasks.

**Please note that no original reference is kept to the parent string/memory, so `WeakRefString` becomes unsafe
once the parent object goes out of scope (i.e. loses a reference to it)**

Internally, a `WeakRefString{T}` holds:

  * `ptr::Ptr{T}`: a pointer to the string data (code unit size is parameterized on `T`)
  * `len::Int`: the number of code units in the string data

See also [`WeakRefStringArray`](@ref)
"""
struct WeakRefString{T} <: AbstractString
    ptr::Ptr{T}
    len::Int # of code units
end

WeakRefString(ptr::Ptr{T}, len) where {T} = WeakRefString(ptr, Int(len))
WeakRefString(t::Tuple{Ptr{T}, Int}) where {T} = WeakRefString(t[1], t[2])

const NULLSTRING = WeakRefString(Ptr{UInt8}(0), 0)
const NULLSTRING16 = WeakRefString(Ptr{UInt16}(0), 0)
const NULLSTRING32 = WeakRefString(Ptr{UInt32}(0), 0)
Base.zero(::Type{WeakRefString{T}}) where {T} = WeakRefString(Ptr{T}(0), 0)

import Base: ==
function ==(x::WeakRefString{T}, y::WeakRefString{T}) where {T}
    x.len == y.len && (x.ptr == y.ptr || ccall(:memcmp, Cint, (Ptr{T}, Ptr{T}, Csize_t),
                                           x.ptr, y.ptr, x.len) == 0)
end
function ==(x::String, y::WeakRefString{T}) where {T}
    sizeof(x) == y.len && ccall(:memcmp, Cint, (Ptr{T}, Ptr{T}, Csize_t),
                                 pointer(x), y.ptr, y.len) == 0
end
==(y::WeakRefString, x::String) = x == y

function Base.hash(s::WeakRefString{T}, h::UInt) where {T}
    h += Base.memhash_seed
    ccall(Base.memhash, UInt, (Ptr{T}, Csize_t, UInt32), s.ptr, s.len, h % UInt32) + h
end

Base.show(io::IO, ::Type{WeakRefString{T}}) where {T} = print(io, "WeakRefString{$T}")
function Base.show(io::IO, x::WeakRefString{T}) where {T}
    print(io, '"')
    print(io, string(x))
    print(io, '"')
    return
end
Base.print(io::IO, s::WeakRefString) = print(io, string(s))
if isdefined(Base, :textwidth)
    Base.textwidth(s::WeakRefString) = textwidth(string(s))
else
    Base.strwidth(s::WeakRefString) = strwidth(string(s))
end

chompnull(x::WeakRefString{T}) where {T} = unsafe_load(x.ptr, x.len) == T(0) ? x.len - 1 : x.len

Base.string(x::WeakRefString) = x == NULLSTRING ? "" : unsafe_string(x.ptr, x.len)
Base.string(x::WeakRefString{UInt16}) = x == NULLSTRING16 ? "" : String(transcode(UInt8, unsafe_wrap(Array, x.ptr, chompnull(x))))
Base.string(x::WeakRefString{UInt32}) = x == NULLSTRING32 ? "" : String(transcode(UInt8, unsafe_wrap(Array, x.ptr, chompnull(x))))

Base.convert(::Type{WeakRefString{UInt8}}, x::String) = WeakRefString(pointer(x), sizeof(x))
Base.convert(::Type{String}, x::WeakRefString) = convert(String, string(x))
Base.String(x::WeakRefString) = string(x)
Base.Symbol(x::WeakRefString{UInt8}) = ccall(:jl_symbol_n, Ref{Symbol}, (Ptr{UInt8}, Int), x.ptr, x.len)

init(::Type{T}, rows) where {T} = fill(zero(T), rows)
if !isdefined(Base, :uninitialized)
    struct Uninitialized end
    const uninitialized = Uninitialized()
    Vector{T}(::Uninitialized, rows) where {T} = Vector{T}(rows)
end
init(::Type{Union{Missing, T}}, rows) where {T} = Vector{Union{Missing, T}}(uninitialized, rows)

"""
A [`WeakRefString`](@ref) container.
Holds the "strong" references to the data pointed by its strings, ensuring that
the referenced memory blocks stay valid during `WeakRefStringArray` lifetime.

Upon indexing an elemnt in a `WeakRefStringArray`, the underlying `WeakRefString` is converted to a proper
Julia `String` type by copying the memory; this ensures safe string processing in the general case. If additional
optimizations are desired, the direct `WeakRefString` elements can be accessed by indexing `A.elements`, where
`A` is a `WeakRefStringArray`.
"""
struct WeakRefStringArray{T<:WeakRefString, N, U} <: AbstractArray{Union{String, U}, N}
    data::Vector{Any}
    elements::Array{Union{T, U}, N}

    WeakRefStringArray(data::Vector{Any}, A::Array{Union{T, Missing}, N}) where {T <: WeakRefString, N} =
        new{T, N, Missing}(data, A)
    WeakRefStringArray(data::Vector{Any}, A::Array{T, N}) where {T <: WeakRefString, N} =
        new{T, N, Union{}}(data, A)
end

WeakRefStringArray(data, ::Type{T}, rows::Integer) where {T} = WeakRefStringArray(Any[data], init(T, rows))
WeakRefStringArray(data, A::Array{T}) where {T <: Union{WeakRefString, Missing}} = WeakRefStringArray(Any[data], A)

wk(A, B::AbstractArray) = WeakRefStringArray(A.data, B)
wk(A, w::WeakRefString) = string(w)
wk(A, ::Missing) = missing

Base.size(A::WeakRefStringArray) = size(A.elements)
Base.getindex(A::WeakRefStringArray, I...) = wk(A, getindex(A.elements, I...))
Base.setindex!(A::WeakRefStringArray{T, N}, v::Missing, i::Int) where {T, N} = setindex!(A.elements, v, i)
Base.setindex!(A::WeakRefStringArray{T, N}, v::Missing, I::Vararg{Int, N}) where {T, N} = setindex!(A.elements, v, I...)
Base.setindex!(A::WeakRefStringArray{T, N}, v::WeakRefString, i::Int) where {T, N} = setindex!(A.elements, v, i)
Base.setindex!(A::WeakRefStringArray{T, N}, v::WeakRefString, I::Vararg{Int, N}) where {T, N} = setindex!(A.elements, v, I...)
Base.setindex!(A::WeakRefStringArray{T, N}, v::String, i::Int) where {T, N} = (push!(A.data, codeunits(v)); setindex!(A.elements, v, i))
Base.setindex!(A::WeakRefStringArray{T, N}, v::String, I::Vararg{Int, N}) where {T, N} = (push!(A.data, codeunits(v)); setindex!(A.elements, v, I...))
Base.resize!(A::WeakRefStringArray, i) = resize!(A.elements, i)

Base.push!(a::WeakRefStringArray{T, 1}, v::Missing) where {T} = (push!(a.elements, v); a)
Base.push!(a::WeakRefStringArray{T, 1}, v::WeakRefString) where {T} = (push!(a.elements, v); a)
function Base.push!(A::WeakRefStringArray{T, 1}, v::String) where T
    push!(A.data, codeunits(v))
    push!(A.elements, v)
    return A
end

function Base.append!(a::WeakRefStringArray{T, 1}, b::WeakRefStringArray{T, 1}) where {T}
    append!(a.data, b.data)
    append!(a.elements, b.elements)
    return a
end

function Base.vcat(a::WeakRefStringArray{T, 1}, b::WeakRefStringArray{T, 1}) where T
    WeakRefStringArray(Any[a.data, b.data], vcat(a.elements, b.elements))
end

end # module


text:
function transform_density!(out::GenericGrayImage, img::GenericGrayImage, edges::AbstractRange, newvals::AbstractVector)
    first_edge, last_edge = first(edges), last(edges)
    first_newval, last_newval = first(newvals), last(newvals)
    inv_step_size = 1/step(edges)
    function transform(val)
        val = gray(val)
        if val >= last_edge
            return last_newval
        elseif val < first_edge
            return first_newval
        else
            index = floor(Int, (val-first_edge)*inv_step_size) + 1
            @inbounds newval = newvals[index]
            return newval
        end
    end
    map!(transform, out, img)
end

function transform_density!(out::GenericGrayImage, img::GenericGrayImage{T}, edges::AbstractRange, newvals::AbstractVector) where T<:Union{N0f8,AbstractGray{N0f8}}
    # When dealing with 8-bit images, we can improve computational performance by precomputing the lookup table
    # for how the intensities transform (there are only 256 calculations of intensities rather than `length(img)`
    # calculations of intensities).
    lookup = Vector{eltype(newvals)}(undef, 256)
    invoke(transform_density!, Tuple{GenericGrayImage,GenericGrayImage,AbstractRange,AbstractVector},
                               lookup, zero(T):eps(T):oneunit(T), edges, newvals)
    map!(out, img) do val
        lookup[uint8val(val)+1]
    end
end

uint8val(x::N0f8) = reinterpret(x)
uint8val(x::AbstractGray) = uint8val(gray(x))

function build_lookup(cdf, minval::T, maxval::T) where T
    first_cdf = first(cdf)
    # Scale the new intensity value to so that it lies in the range [minval, maxval].
    scale = (maxval - minval) / (cdf[end] - first_cdf)
    if T <: Integer
        return T[ceil(minval + (x - first_cdf) * scale) for x in cdf]
    end
    return T[minval + (x - first_cdf) * scale for x in cdf]
end

function construct_pdfs(img::GenericGrayImage, targetimg::AbstractArray, edges::AbstractRange)
    _, histogram = build_histogram(img, edges)
    _, target_histogram = build_histogram(targetimg, edges)
    return edges, histogram / sum(histogram), target_histogram / sum(target_histogram)
end

function construct_pdfs(img::GenericGrayImage, targetimg::AbstractArray, nbins::Integer = 256)
    if eltype(img) <: AbstractGray
        imin, imax = 0, 1
    else
        imin, imax = min(minfinite(img), minfinite(targetimg)), max(maxfinite(img), maxfinite(targetimg))
    end
    edges, histogram = build_histogram(img, nbins, minval = imin, maxval = imax)
    _, target_histogram = build_histogram(targetimg, edges)
    return edges, histogram / sum(histogram), target_histogram / sum(target_histogram)
end

function lookup_icdf(cdf::AbstractArray, targetcdf::AbstractArray)
    lookup_table = zeros(Int, length(cdf))
    i = 1
    for j = 1:length(cdf)
        p = cdf[j]
        while i < length(targetcdf) && targetcdf[i+1] <= p
            i += 1
        end
        lookup_table[j] = i
    end
    lookup_table
end

function cdf2pdf!(pdf::AbstractArray, cdf::AbstractArray)
    pdf[1] = cdf[1]
    for i = length(cdf)-1:-1:2
        pdf[i] = cdf[i] - cdf[i-1]
    end
end


text:
# plot a vector of points
@recipe function f(pts::AbstractVector{<:SVector{N}}) where {N}
    if N == 2
        xx = [pt[1] for pt in pts]
        yy = [pt[2] for pt in pts]
        return xx,yy
    elseif N == 3
        xx = [pt[1] for pt in pts]
        yy = [pt[2] for pt in pts]
        zz = [pt[3] for pt in pts]
        return xx,yy,zz
    else
        notimplemented()
    end
end

@recipe function f(pts::AbstractMatrix{<:SVector})
    vec(pts)
end

# plot a hyperrectangle
@recipe function f(rec::HyperRectangle{N}) where {N}
    seriestype := :path
    label := ""
    if N == 2
        pt1 = rec.low_corner
        pt2 = rec.high_corner
        x1, x2 = pt1[1],pt2[1]
        y1, y2 = pt1[2],pt2[2]
        [x1,x2,x2,x1,x1],[y1,y1,y2,y2,y1]
    elseif N == 3
        seriestype := :path
        pt1 = rec.low_corner
        pt2 = rec.high_corner
        x1, x2 = pt1[1],pt2[1]
        y1, y2 = pt1[2],pt2[2]
        z1, z2 = pt1[3],pt2[3]
        @series [x1,x2,x2,x1,x1],[y1,y1,y2,y2,y1],[z1,z1,z1,z1,z1]
        @series [x1,x2,x2,x1,x1],[y1,y1,y2,y2,y1],[z2,z2,z2,z2,z2]
        @series [x1,x1],[y1,y1],[z1,z2]
        @series [x2,x2],[y1,y1],[z1,z2]
        @series [x2,x2],[y2,y2],[z1,z2]
        @series [x1,x1],[y2,y2],[z1,z2]
    end
end

# plot domain --> plot all of its entities
@recipe function f(Ω::Domain)
    for ent in entities(Ω)
        @series begin
            ent
        end
    end
end

# recipe for parametric line
@recipe function f(ent::ParametricEntity{ReferenceLine};h=0.01)
    par = ent.parametrization
    grid   --> false
    aspect_ratio --> :equal
    s       = 0:h:1
    pts     = [par(v) for v in s]
    x       = [pt[1] for pt in pts]
    y       = [pt[2] for pt in pts]
    x,y
end

@recipe function f(ents::Vector{ParametricEntity{ReferenceLine}};h=0.01)
    grid   --> false
    aspect_ratio --> :equal
    for ent in ents
        @series begin
            ent
        end
    end
end

# recipe for paramatric surface
@recipe function f(ent::ParametricEntity{ReferenceSquare};h=0.1)
    legend --> false
    grid   --> false
    # aspect_ratio --> :equal
    seriestype := :surface
    xrange = 0:h:1
    yrange = 0:h:1
    pts    = [ent.parametrization((x,y)) for x in xrange, y in yrange]
    x      =  [pt[1] for pt in pts]
    y      =  [pt[2] for pt in pts]
    z      =  [pt[3] for pt in pts]
    x,y,z
end

# recipe for parametric body
@recipe function f(bdy::AbstractParametricBody)
    label --> ""
    # aspect_ratio --> :equal
    for patch in boundary(bdy)
        @series begin
            patch
        end
    end
end

# recipe for many parametric bodies
@recipe function f(bdies::Vector{<:AbstractParametricBody})
    aspect_ratio --> :equal
    for bdy in bdies
        @series begin
            bdy
        end
    end
end

@recipe function f(mesh::NystromMesh)
    label --> ""
    grid   --> false
    aspect_ratio --> :equal
    for (E,els) in mesh.elements
        for el in els
            @series begin
                el
            end
        end
    end
end

@recipe function f(mesh::GenericMesh,Ω::Domain)
    view(mesh,Ω)
end

@recipe function f(mesh::SubMesh)
    label --> ""
    grid   --> false
    aspect_ratio --> :equal
    for E in etypes(mesh)
        for el in ElementIterator(mesh,E)
            @series begin
                el
            end
        end
    end
end

# FIXME: write better recipes
@recipe function f(el::LagrangeTriangle)
    label --> ""
    vtx = nodes(el)
    for n in 1:3
        is = n
        ie = 1 + (n%3)
        @series begin
            [vtx[is],vtx[ie]]
        end
    end
end

@recipe function f(el::LagrangeRectangle)
    label --> ""
    vtx = nodes(el)
    for n in 1:4
        is = n
        ie = 1 + (n%4)
        @series begin
            [vtx[is],vtx[ie]]
        end
    end
end

@recipe function f(el::LagrangeLine)
    vtx = nodes(el)
    [vtx[1],vtx[2]]
end

@recipe function f(el::AbstractElement{ReferenceLine};npts=10)
    label --> ""
    grid   --> false
    s       =  LinRange(0,1,npts)
    pts     = [el(v) for v in s]
    x       = [pt[1] for pt in pts]
    y       = [pt[2] for pt in pts]
    return x,y
end

@recipe function f(els::NTuple{<:Any,<:AbstractElement{ReferenceLine}};npts=10)
    label --> ""
    grid   --> false
    for el in els
        @series begin
            el
        end
    end
end

@recipe function f(ent::AbstractElement{ReferenceSquare};h=0.1)
    legend --> false
    grid   --> true
    seriesalpha --> 0.1
    # aspect_ratio --> :equal
    seriestype --> :wireframe
    xrange = 0:h:1
    yrange = 0:h:1
    pts    = [ent((x,y)) for x in xrange, y in yrange]
    x      =  [pt[1] for pt in pts]
    y      =  [pt[2] for pt in pts]
    z      =  [pt[3] for pt in pts]
    x,y,z
end

max_stars_repo_name: WaveProp/WaveProp

max_issues_repo_name: WaveProp/WaveProp

max_forks_repo_name: WaveProp/WaveProp


text:
mutable struct DLNode{T}
    value::T
    pred::Union{DLNode{T}, Nothing}
    succ::Union{DLNode{T}, Nothing}
    DLNode(v) = new{typeof(v)}(v, nothing, nothing)
end

function insertpost(prevnode, node)
    succ = prevnode.succ
    prevnode.succ = node
    node.pred = prevnode
    node.succ = succ
    if succ != nothing
        succ.pred =  node
    end
    node
end

function insertpre(postnode, node)
    pred = postnode.pred
    postnode.pred = node
    node.succ = postnode
    node.pred = pred
    if pred != nothing
        pred.succ = node
    end
    node
end

function delete(nd)
    if nd != nothing
        pred = nd.pred
        succ = nd.succ
        if pred != nothing pred.succ = succ end
        if succ != nothing succ.pred = pred end
    end
    nothing
end

first(nd) = (while nd.pred != nothing nd = nd.prev end; nd)
last(nd) = (while nd.succ != nothing nd = nd.succ end; nd)

function printconnected(nd; fromtail = false)
    if fromtail
        nd = last(nd)
        print(nd.value)
        while nd.pred != nothing
            nd = nd.pred
            print(" -> $(nd.value)")
        end
    else
        nd = first(nd)
        print(nd.value)
        while nd.succ != nothing
            nd = nd.succ
            print(" -> $(nd.value)")
        end
    end
    println()
end

node1 = DLNode(1)
node2 = DLNode(2)
node3 = DLNode(3)
node4 = DLNode(4)
insertpost(node1, node2)
insertpre(node2, node4)
insertpost(node2, node3)
println("First value is ", first(node1).value, " and last value is ", last(node1).value)
print("From beginning to end: "); printconnected(node1)
delete(node4)
print("From end to beginning post deletion: "); printconnected(node1, fromtail = true)

max_stars_repo_name: ethansaxenian/RosettaDecode

max_issues_repo_name: ethansaxenian/RosettaDecode

max_forks_repo_name: ethansaxenian/RosettaDecode


text:
# coding: utf-8

# In[1]:

#this is a test on how recursion with memory reduces time complexity

# In[2]:

#this is the fibonacci recursion function without memory
#it is basically algorithm 101 for any coding language
function fib(n)

    if n==1
        return 1
    elseif n==2
        return 1
    elseif n<=0
        printstyled("Invalid input",color=:red)
        return
    else
        return fib(n-1)+fib(n-2)

    end
end

# In[3]:

#i need a global dictionary to do the memorization
#or i can change function mmz(n) into mmz(n,memoization)
global memoization=Dict(1=>1,2=>1)

# In[4]:

#mmz(n) is recursion with memory
#everytime we do the calculation, we store it in the dictionary
#i denote the key as the n th fibonacci number
#the value as the number itself
#if we can find the key in dictionary
#we simply return the value
#if not, we compute and update the dictionary then return the value
function mmz(n)

    if n<=0
        printstyled("Invalid input",color=:red)
        return
    end

    if !(n in keys(memoization))
        global memoization[n]=mmz(n-1)+mmz(n-2)
    end

    return memoization[n]
end

# In[5]:

#using ijulia inline magic
#equivalent to %timeit in ipython
#0.000093 seconds
@time begin
    fib(20)
end

# In[6]:

#0.019746 seconds (5.60 k allocations: 335.035 KiB)
#it seems to be slower
#but if we compute mmz(30) now
#its much faster than fib(30)
@time begin
    mmz(20)
end

# In[7]:

#0.021588 seconds
@time begin
    fib(30)
end

# In[8]:

#0.000023 seconds (52 allocations: 832 bytes)
@time begin
    mmz(30)
end

max_stars_repo_path: fibonacci with memoization.jl

max_issues_repo_path: fibonacci with memoization.jl

max_forks_repo_path: fibonacci with memoization.jl


text:
import LinearAlgebra: \, ldiv!, mul!, norm, LU, factorize, Factorization, dot

using LinearMaps

abstract type AbstractPotentialFlowSystem{TU} end

"""
$(TYPEDEF)

Defines an immersed boundary Poisson system that solves for the streamfunction and bound vortex sheet strength, given a streamfunction boundary condition on immersed boundary points.

# Fields

$(TYPEDFIELDS)
"""
struct IBPoisson{TU,TF} <: AbstractPotentialFlowSystem{TU}
    """L: Discrete Laplacian."""
    L::CartesianGrids.Laplacian
    """R: Regularization operator."""
    R::RegularizationMatrix{TU,TF}
    """E: Interpolation operator."""
    E::InterpolationMatrix{TU,TF}
    """Sfact: factorized version of the Schur complement."""
    Sfact::LU # Factorization

    """Buffers"""
    _A⁻¹r₁::TU
    _B₁ᵀf::TU
    _f_buf::TF
end

function IBPoisson(L::CartesianGrids.Laplacian, R::RegularizationMatrix{TU,TF}, E::InterpolationMatrix{TU,TF}) where {TU,TF}
    _A⁻¹r₁ = TU()
    _B₁ᵀf = TU()
    _f_buf = TF()

    Rmap = LinearMap(x->vec(R*TF(x)),length(_A⁻¹r₁),length(_f_buf));
    Emap = LinearMap(x->vec(E*TU(x)),length(_f_buf),length(_A⁻¹r₁))
    L⁻¹map = LinearMap(x->vec(L\TU(x)),length(_A⁻¹r₁));
    Smap = -Emap*L⁻¹map*Rmap
    Sfact = factorize(Matrix(Smap))

    IBPoisson{TU,TF}(L, R, E, Sfact, _A⁻¹r₁, _B₁ᵀf, _f_buf)
end

"""
$(TYPEDEF)

Defines an immersed boundary Poisson system with a single type of constraints on the bound vortex sheet strength.

# Fields

$(TYPEDFIELDS)
"""
struct ConstrainedIBPoisson{Nb,TU,TF} <: AbstractPotentialFlowSystem{TU}
    """ibp: Immersed boundary poisson inner system that solves for the streamfunction and bound vortex sheet strength, given a streamfunction boundary condition on immersed boundary points."""
    ibp::IBPoisson{TU,TF}
    """B₁ᵀ₂cols: Columns of the second partition of the B₁ᵀ matrix for the outer system."""
    B₁ᵀ₂cols::Vector{TF}
    """B₂₂rows: Rows of the second partition of the B₂ matrix for the outer system."""
    B₂₂rows::Vector{TF}
    """Soutfact: Factorized version of the Schur complement for the outer system."""
    Soutfact::Union{Factorization,Float64}

    """Buffers"""
    _f_buf::TF
    _ψ₀_buf::Vector{Float64}
    _sol_buf::IBPoissonSolution{TU,TF}
end

function ConstrainedIBPoisson(L::CartesianGrids.Laplacian, R::RegularizationMatrix{TU,TF}, E::InterpolationMatrix{TU,TF}, B₁ᵀ₂cols::Vector{TF}, B₂₂rows::Vector{TF}) where {TU,TF}
    _f_buf = TF()
    Nb = length(B₁ᵀ₂cols) # number of bodies
    ibp = IBPoisson(L,R,E)
    S₀ = Matrix{Float64}(undef, Nb, Nb)
    _S⁻¹B₁ᵀ₂ = zeros(length(_f_buf), Nb)
    Soutfact = _computeSoutfact(S₀, ibp.Sfact, B₁ᵀ₂cols, B₂₂rows, _S⁻¹B₁ᵀ₂)
    _ψ₀_buf = zeros(Nb)
    _sol_buf = IBPoissonSolution{TU,TF}(TU(),TF())
    ConstrainedIBPoisson{Nb,TU,TF}(ibp, B₁ᵀ₂cols, B₂₂rows, Soutfact, _f_buf, _ψ₀_buf, _sol_buf)
end

"""
$(TYPEDEF)

Defines an immersed boundary Poisson system with constraints on the bound vortex sheet strength that regularize up to one edge per body.

# Fields

$(TYPEDFIELDS)
"""
struct SteadyRegularizedIBPoisson{Nb,Ne,TU,TF} <: AbstractPotentialFlowSystem{TU}
    """cibp: Constrained immersed boundary Poisson system."""
    cibp::ConstrainedIBPoisson{Nb,TU,TF}
    """one_vec: Array of vectors with the i-th entry of the j-th vector equal to one if i is an index of the points that belong to the j-th body and zero otherwise."""
    one_vec::Vector{TF}
    """e_vec: Array of unit vectors with the inner product of the k-th unit vector and vortex sheet strength equal to the value of the vortex sheet strength at the k-th edge."""
    e_vec::Vector{TF}
    """f₀: Unregularized bound vortex sheet strength due to a uniform unit boundary condition for the streamfunction on all immersed bodies."""
    f₀::TF
    """f₀_vec: Array of vectors with the i-th entry of the j-th vector equal to f₀[i] if i is an index of the points that belong to the j-th body and zero otherwise."""
    f₀_vec::Vector{TF}
end

function SteadyRegularizedIBPoisson(L::CartesianGrids.Laplacian, R::RegularizationMatrix{TU,TF}, E::InterpolationMatrix{TU,TF}, one_vec::Vector{TF}, e_vec::Vector{TF}) where {TU,TF}

    Nb = length(one_vec) # number of bodies
    Ne = length(e_vec) # number of edges

    _TU_buf = TU()

    f₀ = TF()
    _TF_buf = TF()
    _TF_buf .= 1.0
    ibp = IBPoisson(L, R, E)
    ldiv!(IBPoissonSolution(_TU_buf,f₀),ibp,IBPoissonRHS(_TU_buf,_TF_buf), onlyf=true, zerow=true)
    f₀_vec = [f₀.*one_vec[i] for i=1:Nb]

    R̃ = deepcopy(R)
    R̃.M .= R̃.M*Diagonal(f₀)
    cibp = ConstrainedIBPoisson(L, R̃, E, one_vec, e_vec)
    SteadyRegularizedIBPoisson{Nb,Ne,TU,TF}(cibp, one_vec, e_vec, f₀, f₀_vec)
end

"""
$(TYPEDEF)

Defines an immersed boundary Poisson system that, in addition to solving for the streamfunction and bound vortex sheet strength, solves for N point vortex strengths to satisfy N constraints on the bound vortex sheet strength that regularize specified edges.

# Fields

$(TYPEDFIELDS)
"""
struct UnsteadyRegularizedIBPoisson{Nb,Ne,TU,TF} <: AbstractPotentialFlowSystem{TU}
    """ibp: Immersed boundary poisson inner system with re-scaled regularization operator that solves for the streamfunction and bound vortex sheet strength, given a streamfunction boundary condition on immersed boundary points."""
    ibp::IBPoisson{TU,TF}
    """one_vec: Array of vectors with the i-th entry of the j-th vector equal to one if i is an index of the points that belong to the j-th body and zero otherwise."""
    one_vec::Vector{TF}
    """e_vec: Array of unit vectors with the inner product of the i-th unit vector and vortex sheet strength equal to the value of the vortex sheet strength at the i-th edge."""
    e_vec::Vector{TF}
    """vidx_vec: Array of vectors with the i-th vector containing the indices of the vortices that are shedded from the i-th body."""
    vidx_vec::Vector{Vector{Int}}
    """f₀: Unregularized vortex sheet strength due to a uniform unit boundary condition for the streamfunction on all immersed bodies."""
    f₀::TF
    """f₀_vec: Array of vectors with the i-th entry of the j-th vector equal to f₀[i] if i is an index of the points that belong to the j-th body and zero otherwise."""
    f₀_vec::Vector{TF}
    """f̃₀_vec: Array of unregularized bound vortex sheet strengths, obtained using the re-scaled regularization operator, with the i-th entry the bound vortex sheet strength due to a uniform unit boundary condition for the streamfunction on the i-th immersed body and zero everywhere else."""
    f̃₀_vec::Vector{TF}
    """f̃_vec: Array of unregularized bound vortex sheet strengths, obtained using the re-scaled regularization operator, with the i-th entry the bound vortex sheet strength due to the i-th vorticity field in `d_vec` and a zero boundary condition for the streamfunction on all bodies."""
    f̃_vec::Vector{TF}
    """Γ₀: The circulation of the of the bound vortex sheet strength f₀."""
    Γ₀::Float64
    """d_vec: Array of vorticity fields with the i-th entry the vorticity field due to the i-th point vortex that is used for regularizing the edges with its strength set to one and the other point vortex strengths set to zero."""
    d_vec::Vector{TU}
    """activef̃lim_vec: Array of f̃ values that are used as the right-hand side for the constraints on f̃ when solving the system. These values are chosen by comparing the unconstrained f̃ from the inner system to the limits provided in `rhs` to `ldiv!` and selecting the limit that is exceeded or setting it to `Inf` if no limit is exceeded."""
    activef̃lim_vec::Vector{f̃Limit}
    """Sout: Schur complement for the outer system."""
    Sout::Matrix{Float64}

    """Buffers"""
    _TU_buf::TU
    _f_buf::TF
    _r₂_buf::Vector{Float64}
    _y_buf::Vector{Float64}
end

function UnsteadyRegularizedIBPoisson(L::CartesianGrids.Laplacian, R::RegularizationMatrix{TU,TF}, E::InterpolationMatrix{TU,TF}, one_vec::Vector{TF}, e_vec::Vector{TF}, vidx_vec::Vector{Vector{Int}}) where {TU,TF}
    Nb = length(one_vec) # number of bodies
    Ne = length(e_vec) # number of edges

    ibp = IBPoisson(L, R, E)

    _TU_buf = TU()
    _f_buf = TF()

    f₀ = TF()

    _f_buf .= 1.0
    ldiv!(IBPoissonSolution(_TU_buf,f₀),ibp,IBPoissonRHS(_TU_buf,_f_buf), onlyf=true, zerow=true)

    f₀_vec = [f₀.*one_vec[i] for i=1:Nb]

    R̃ = deepcopy(R)
    R̃.M .= R̃.M*Diagonal(f₀)
    ibp = IBPoisson(L, R̃, E)

    f̃₀_vec = [TF() for i=1:Nb]
    _computef₀_vec!(f̃₀_vec, ibp, one_vec, _TU_buf)

    d_vec = [TU() for i=1:Ne]
    f̃_vec = [TF() for i=1:Ne]
    Γ₀ = sum(f₀)
    activef̃lim_vec = zeros(Ne)
    Sout = zeros(Nb+Ne,Nb+Ne)
    _r₂_buf = zeros(Nb+Ne)
    _y_buf = zeros(Nb+Ne)

    UnsteadyRegularizedIBPoisson{Nb,Ne,TU,TF}(ibp, one_vec, e_vec, vidx_vec, f₀, f₀_vec, f̃₀_vec, f̃_vec, Γ₀, d_vec, activef̃lim_vec, Sout, _TU_buf, _f_buf, _r₂_buf, _y_buf)
end

function ldiv!(sol::TS, sys::IBPoisson{TU,TF}, rhs::TR; zerow=false, zeroψb=false, useprovidedf=false, onlyf=false) where {TU,TF,TS<:AbstractIBPoissonSolution,TR<:AbstractIBPoissonRHS}

    if useprovidedf # in this case we use the f provided in sol.f and only calculate ψ = L⁻¹(-w-Rf)
        mul!(sys._B₁ᵀf, sys.R, sol.f)
        sys._B₁ᵀf .= .-rhs.w .- sys._B₁ᵀf
        ldiv!(sol.ψ, sys.L, sys._B₁ᵀf)
    else
        if zerow
            sys._f_buf .= 0.0
        else
            ldiv!(sys._A⁻¹r₁, sys.L, rhs.w) # -A⁻¹r₁ (note minus sign because w and not -w is used)
            mul!(sys._f_buf, sys.E, sys._A⁻¹r₁) # -B₂A⁻¹r₁ (note minus sign)
        end

        if !zeroψb
            sys._f_buf .= rhs.ψb .+ sys._f_buf # r₂ - B₂A⁻¹r₁
        end

        ldiv!(sol.f, sys.Sfact, sys._f_buf) # S⁻¹(r₂ - B₂A⁻¹r₁)

        if !onlyf # ψ = L⁻¹(-w-Rf)
            mul!(sys._B₁ᵀf, sys.R, sol.f)
            ldiv!(sol.ψ, sys.L, sys._B₁ᵀf)
            sol.ψ .= sys._A⁻¹r₁ .- sol.ψ # L⁻¹(-w) was already calculated earlier and stored in sys._A⁻¹r₁
        end
    end

    return sol
end

function ldiv!(sol::ConstrainedIBPoissonSolution{TU,TF}, sys::ConstrainedIBPoisson{Nb,TU,TF}, rhs::ConstrainedIBPoissonRHS) where {Nb,TU,TF,TS}

    # fstar = S⁻¹(ψb+EL⁻¹w)
    ldiv!(sol, sys.ibp, rhs, onlyf=true)

    # ψ₀ = S₀\(Γb-B₂₂*fstar)
    for b in 1:Nb
        sys._ψ₀_buf[b] = dot(sys.B₂₂rows[b]',sol.f)
    end
    sys._ψ₀_buf .= rhs.fconstraintRHS .- sys._ψ₀_buf
    ldiv!(sol.ψ₀, sys.Soutfact, sys._ψ₀_buf)

    # f = fstar - S⁻¹(B₁ᵀ₂*ψ₀)
    sys._f_buf .= 0.0
    _addlincomboto!(sys._f_buf, sol.ψ₀, sys.B₁ᵀ₂cols)
    ldiv!(sys.ibp.Sfact, sys._f_buf)
    sol.f .-= sys._f_buf

    # ψ = L⁻¹(-w - Rf)
    ldiv!(sol, sys.ibp, rhs, useprovidedf=true)

end

function ldiv!(sol::ConstrainedIBPoissonSolution{TU,TF}, sys::SteadyRegularizedIBPoisson{Nb,Ne,TU,TF}, rhs::ConstrainedIBPoissonRHS) where {Nb,Ne,TU,TF}

    ldiv!(sol, sys.cibp, rhs)
    sol.f .*= sys.f₀

end

function ldiv!(sol::ConstrainedIBPoissonSolution{TU,TF}, sys::UnsteadyRegularizedIBPoisson{Nb,Ne,TU,TF}, rhs::UnsteadyRegularizedIBPoissonRHS{TU,TF}) where {Nb,Ne,TU,TF}

    # f̃star = S̃⁻¹(ψb+EL⁻¹w)
    ldiv!(sol, sys.ibp, rhs, onlyf=true)

    for k in 1:Ne
        # eq 2.47
        ldiv!(IBPoissonSolution(sys._TU_buf, sys.f̃_vec[k]), sys.ibp, IBPoissonRHS(sys.d_vec[k],sys._f_buf), zeroψb=true, onlyf=true)
        # eq 2.60
        sys.activef̃lim_vec[k] = _findactivef̃limit(sys.e_vec[k], sol.f, rhs.f̃lim_vec[k])
    end

    # The shedding edges are the ones for which the active f̃ limit is not Inf
    sheddingedges = [e for e in 1:length(sys.activef̃lim_vec) if sys.activef̃lim_vec[e] != Inf]

    # eq 2.49, eq 2.56, eq 2.61
    _computeδΓandψ₀!(sol,sys,rhs,sheddingedges)

    # f̃ = f̃star + Σf̃ᵢδΓᵢ - Σf̃₀ⱼψ₀ⱼ
    _addlincomboto!(sol.f, sol.δΓ_vec, sys.f̃_vec)
    sys._f_buf .= 0.0
    _addlincomboto!(sys._f_buf, sol.ψ₀, sys.f̃₀_vec)
    sol.f .-= sys._f_buf

    # ψ = L⁻¹(-w - ΣdᵢδΓᵢ - Rf)
    sys._TU_buf .= rhs.w
    _addlincomboto!(sys._TU_buf, sol.δΓ_vec, sys.d_vec)
    ldiv!(sol, sys.ibp, IBPoissonRHS(sys._TU_buf, sys._f_buf), useprovidedf=true)

    sol.f .*= sys.f₀

end

function _addlincomboto!(s,c,A)
    for i in 1:length(c)
        s .+= c[i].*A[i]
    end
end

function _computeSoutfact(S₀,Sfact,B₁ᵀ₂cols,B₂₁rows,_S⁻¹B₁ᵀ₂)
    for c in 1:length(B₁ᵀ₂cols)
        ldiv!(view(_S⁻¹B₁ᵀ₂,:,c),Sfact,B₁ᵀ₂cols[c])
    end

    for r in 1:length(B₂₁rows)
        mul!(view(S₀,r:r,:),B₂₁rows[r]',_S⁻¹B₁ᵀ₂) # need to find better solution for this indexing
    end

    S₀ .= .-S₀

    Soutfact = factorize(S₀)
end

function ldiv!(x::ScalarData ,A::LU, f::ScalarData)
    ldiv!(x.data,A,f.data)
end

function ldiv!(A::LU, f::ScalarData)
    ldiv!(A,f.data)
end

"""
$(SIGNATURES)

Computes the vortex strengths for the edges `sheddingedges` and ψ₀ values and stores them in the `δΓ_vec` and `ψ₀` fields of `sol`. If there are edges that are not in `sheddingedges`, their corresponding vortex strengths are set to zero.
"""
function _computeδΓandψ₀!(sol::ConstrainedIBPoissonSolution{TU,TF}, sys::UnsteadyRegularizedIBPoisson{Nb,Ne,TU,TF}, rhs::UnsteadyRegularizedIBPoissonRHS{TU,TF}, sheddingedges::Vector{Int}) where {Nb,Ne,TU,TF}
    sys.Sout .= 0.0
    for i in 1:Ne
        if i in sheddingedges
            for b in 1:Nb
                sys.Sout[i,b] = dot(sys.e_vec[i], sys.f̃₀_vec[b])
            end
            for e in 1:Ne
                sys.Sout[i,end-Ne+e] = -dot(sys.e_vec[i], sys.f̃_vec[e])
            end
            sys._r₂_buf[i] = sys.activef̃lim_vec[i] - dot(sys.e_vec[i], sol.f)
        else
            sys.Sout[i,end-Ne+i] = 1.0
            sys._r₂_buf[i] = 0.0
        end
    end
    for i in 1:Nb
        for b in 1:Nb
            sys.Sout[end-Nb+i,b] = dot(sys.f₀_vec[i], sys.f̃₀_vec[b])
        end
        for e in 1:Ne
            sys.Sout[end-Nb+i,end-Ne+e] = -dot(sys.f₀_vec[i], sys.f̃_vec[e])
        end
        sys._r₂_buf[end-Nb+i] = -rhs.Γw_vec[i] - dot(sys.f₀_vec[i], sol.f)
    end
    sys.Sout .= .-sys.Sout
    for i in 1:Nb
        for vid in sys.vidx_vec[i]
            sys.Sout[end-Nb+i,end-Ne+vid] += 1.0
        end
    end

    ldiv!(sys._y_buf, factorize(sys.Sout), sys._r₂_buf)
    sol.ψ₀ .= sys._y_buf[1:Nb]
    sol.δΓ_vec .= sys._y_buf[Nb+1:end]
end

function _findactivef̃limit(e::AbstractVector, f̃::TF, f̃lim::f̃Limits) where {TF}

    if e'*f̃ < f̃lim.min
        activef̃lim = f̃lim.min
    elseif e'*f̃ > f̃lim.max
        activef̃lim = f̃lim.max
    else
        activef̃lim = Inf
    end

    return activef̃lim
end

function _computef₀_vec!(f₀_vec::Vector{TF}, ibp::IBPoisson{TU,TF}, one_vec::Vector{TF}, TU_buf::TU) where {TU,TF}
    for b in 1:length(f₀_vec)
        ldiv!(IBPoissonSolution(TU_buf,f₀_vec[b]),ibp,IBPoissonRHS(TU_buf,one_vec[b]), onlyf=true, zerow=true)
    end
end


text:
@propagate_inbounds function ArrayInterface.to_index(::IndexAxis, axis, arg::AbstractRange{I}) where {I<:Integer}
    @boundscheck if !checkindex(Bool, axis, arg)
        throw(BoundsError(axis, arg))
    end
    return arg
end
@propagate_inbounds function ArrayInterface.to_index(::IndexAxis, axis, arg::Integer)
    @boundscheck if !checkindex(Bool, axis, arg)
        throw(BoundsError(axis, arg))
    end
    return arg
end
@propagate_inbounds function ArrayInterface.to_index(::IndexAxis, axis, arg::AbstractArray{Bool})
    @boundscheck if !checkindex(Bool, axis, arg)
        throw(BoundsError(axis, arg))
    end
    return eachindex(axis)[arg]
end
@propagate_inbounds function ArrayInterface.to_index(::IndexAxis, axis, arg::AbstractArray{I}) where {I<:Integer}
    @boundscheck if !checkindex(Bool, axis, arg)
        throw(BoundsError(axis, arg))
    end
    return arg
end
ArrayInterface.to_index(::IndexAxis, axis, ::Colon) = indices(axis)

#= TODO delete these
@propagate_inbounds function ArrayInterface.to_index(::IndexAxis, axis, arg::AbstractUnitRange{I}) where I<:Integer
    return to_index(parent(axis), arg)
end
@propagate_inbounds function _to_index(axis, arg::CartesianIndex)
    @boundscheck checkbounds(axis, arg)
    return arg
end
=#
@propagate_inbounds function ArrayInterface.to_index(::IndexAxis, axis, arg::Union{<:Equal,Approx})
    idx = findfirst(arg, keys(axis))
    @boundscheck if idx isa Nothing
        throw(BoundsError(axis, arg))
    end
    return Int(@inbounds(eachindex(axis)[idx]))
end
@propagate_inbounds function ArrayInterface.to_index(::IndexAxis, axis, arg)
    ks = keys(axis)
    if arg isa keytype(axis)
        idx = find_first(==(arg), ks)
    else
        idx = find_first(==(keytype(axis)(arg)), ks)
    end
    @boundscheck if idx isa Nothing
        throw(BoundsError(axis, (arg,)))
    end
    # if firstindex of kas is not the same as first of parent(axis)
    p = parent(axis)
    kindex = firstindex(ks)
    pindex = first(p)
    return Int(@inbounds(p[idx + (pindex - kindex)]))
end

@propagate_inbounds function ArrayInterface.to_index(::IndexAxis, axis, arg::Interval)
    idx = find_all(in(arg), keys(axis))
    return @inbounds(indices(axis)[idx])  # FIXME
end

# if it's not an `Axis` then there aren't non-integer keys
@propagate_inbounds function ArrayInterface.to_index(::IndexAxis, axis, arg::Function)
    return @inbounds(eachindex(axis)[find_all(arg, keys(axis))])
end

@propagate_inbounds function ArrayInterface.to_index(::IndexAxis, axis, arg::AbstractArray)
    return map(arg_i -> to_index(axis, arg_i), arg)
end

function ArrayInterface.to_index(::IndexAxis, axis, arg::AbstractRange)
    if typeof(arg) <: keytype(axis)
        inds = find_all_in(arg, keys(axis))
        # if `inds` is same length as `arg` then all of `arg` was found and is inbounds
        @boundscheck if length(inds) != length(arg)
            throw(BoundsError(axis, arg))
        end
        return @inbounds(eachindex(axis)[idx])
    else
        # FIXME SubAxis will have problems with this if the parent axis keys selected extend
        # beyond the
        return to_index(parent(axis), arg)
    end
end


text:
function statements(ci::CreateInfo)
    stype_arg = stypes[ci.sdef.name]
    args = collect(keys(ci.sdef.fields))
    args_replaced = args |> Map(x -> replace_count_argument(args, x)) |> collect
    [Statement("create_info = $(ci.sdef.name)($stype_arg, $(join_args(args_replaced)))", :create_info, args_replaced),
    Statement("create_info_ptr = Ref(create_info)", :create_info_ptr, [:create_info])]
end

function statements(create::Create)
    type = create.type # e.g. Device
    ptr_type_vk = "p_ " * convert(SnakeCaseLower, "$type ").value # e.g. p_device
    ptr_type = convert(SnakeCaseLower, "$type ").value * "_ptr" # e.g. device_ptr
    s1 = Statement("$ptr_type_sym = Ref{$type }()", ptr_type, [])
    args = argnames(create.f.signature)
    args_with_ci = args |> Map(x -> convert(SnakeCaseLower, CamelCaseLower("$x")).value) |> Map(x -> startswith(x, "p_") ? x * "_ptr" : x) |> collect
    s2 = Statement("@check $(create.f.name)($(join_args(args_with_ci)))", nothing, ["create_info", ptr_type])
    [s1, s2]
end

function statements(f::Finalizer)
    args = argnames(f.f.signature)
    @assert f.x ∈ args "Finalized variable $(f.x) not in method arguments"
    finalizer_args = args |> Map(x -> x == f.x ? "x" : x) |> collect
    [Statement("Base.finalizer(x -> $(f.f.name)($(join_args(finalizer_args))), $(f.x))", nothing, filter(x -> x != f.x, finalizer_args))]
end

function statements(pats::AbstractArray{<: Pattern})
    stmts = Statement[]
    for pat ∈ pats
        append!(stmts, statements(pat))
    end
end


text:
#=
This file is auto-generated. Do not edit.
=#

#! format: off

"""
    mutable struct STAB1 <: PSS
        KT::Float64
        T::Float64
        T1T3::Float64
        T3::Float64
        T2T4::Float64
        T4::Float64
        H_lim::Float64
        ext::Dict{String, Any}
        states::Vector{Symbol}
        n_states::Int
        states_types::Vector{StateTypes}
        internal::InfrastructureSystemsInternal
    end

Speed-Sensitive Stabilizing Model

# Arguments
- `KT::Float64`: K/T for washout filter, validation range: `(0, nothing)`, action if invalid: `warn`
- `T::Float64`: Time constant for washout filter, validation range: `(0.01, nothing)`, action if invalid: `warn`
- `T1T3::Float64`: Time constant division T1/T3, validation range: `(0, nothing)`
- `T3::Float64`: Time constant, validation range: `(0.01, nothing)`, action if invalid: `warn`
- `T2T4::Float64`: Time constant division T2/T4, validation range: `(0, nothing)`, action if invalid: `warn`
- `T4::Float64`: Time constant, validation range: `(0.01, nothing)`, action if invalid: `warn`
- `H_lim::Float64`: PSS output limit, validation range: `(0, 0.5)`, action if invalid: `warn`
- `ext::Dict{String, Any}`
- `states::Vector{Symbol}`: The states are:
	x_p1: washout filter,
	x_p2: T1/T2 lead-lag block,
	x_p3: T3/T4 lead-lag block,
- `n_states::Int`: STAB1 has 3 states
- `states_types::Vector{StateTypes}`: STAB1 has 3 differential states
- `internal::InfrastructureSystemsInternal`: power system internal reference, do not modify
"""
mutable struct STAB1 <: PSS
    "K/T for washout filter"
    KT::Float64
    "Time constant for washout filter"
    T::Float64
    "Time constant division T1/T3"
    T1T3::Float64
    "Time constant"
    T3::Float64
    "Time constant division T2/T4"
    T2T4::Float64
    "Time constant"
    T4::Float64
    "PSS output limit"
    H_lim::Float64
    ext::Dict{String, Any}
    "The states are:
	x_p1: washout filter,
	x_p2: T1/T2 lead-lag block,
	x_p3: T3/T4 lead-lag block,"
    states::Vector{Symbol}
    "STAB1 has 3 states"
    n_states::Int
    "STAB1 has 3 differential states"
    states_types::Vector{StateTypes}
    "power system internal reference, do not modify"
    internal::InfrastructureSystemsInternal
end

function STAB1(KT, T, T1T3, T3, T2T4, T4, H_lim, ext=Dict{String, Any}(), )
    STAB1(KT, T, T1T3, T3, T2T4, T4, H_lim, ext, [:x_p1, :x_p2, :x_p3], 3, [StateTypes.Differential, StateTypes.Differential, StateTypes.Differential], InfrastructureSystemsInternal(), )
end

function STAB1(; KT, T, T1T3, T3, T2T4, T4, H_lim, ext=Dict{String, Any}(), states=[:x_p1, :x_p2, :x_p3], n_states=3, states_types=[StateTypes.Differential, StateTypes.Differential, StateTypes.Differential], internal=InfrastructureSystemsInternal(), )
    STAB1(KT, T, T1T3, T3, T2T4, T4, H_lim, ext, states, n_states, states_types, internal, )
end

# Constructor for demo purposes; non-functional.
function STAB1(::Nothing)
    STAB1(;
        KT=0,
        T=0.01,
        T1T3=0,
        T3=0.01,
        T2T4=0,
        T4=0.01,
        H_lim=0,
        ext=Dict{String, Any}(),
    )
end

"""Get [`STAB1`](@ref) `KT`."""
get_KT(value::STAB1) = value.KT
"""Get [`STAB1`](@ref) `T`."""
get_T(value::STAB1) = value.T
"""Get [`STAB1`](@ref) `T1T3`."""
get_T1T3(value::STAB1) = value.T1T3
"""Get [`STAB1`](@ref) `T3`."""
get_T3(value::STAB1) = value.T3
"""Get [`STAB1`](@ref) `T2T4`."""
get_T2T4(value::STAB1) = value.T2T4
"""Get [`STAB1`](@ref) `T4`."""
get_T4(value::STAB1) = value.T4
"""Get [`STAB1`](@ref) `H_lim`."""
get_H_lim(value::STAB1) = value.H_lim
"""Get [`STAB1`](@ref) `ext`."""
get_ext(value::STAB1) = value.ext
"""Get [`STAB1`](@ref) `states`."""
get_states(value::STAB1) = value.states
"""Get [`STAB1`](@ref) `n_states`."""
get_n_states(value::STAB1) = value.n_states
"""Get [`STAB1`](@ref) `states_types`."""
get_states_types(value::STAB1) = value.states_types
"""Get [`STAB1`](@ref) `internal`."""
get_internal(value::STAB1) = value.internal

"""Set [`STAB1`](@ref) `KT`."""
set_KT!(value::STAB1, val) = value.KT = val
"""Set [`STAB1`](@ref) `T`."""
set_T!(value::STAB1, val) = value.T = val
"""Set [`STAB1`](@ref) `T1T3`."""
set_T1T3!(value::STAB1, val) = value.T1T3 = val
"""Set [`STAB1`](@ref) `T3`."""
set_T3!(value::STAB1, val) = value.T3 = val
"""Set [`STAB1`](@ref) `T2T4`."""
set_T2T4!(value::STAB1, val) = value.T2T4 = val
"""Set [`STAB1`](@ref) `T4`."""
set_T4!(value::STAB1, val) = value.T4 = val
"""Set [`STAB1`](@ref) `H_lim`."""
set_H_lim!(value::STAB1, val) = value.H_lim = val
"""Set [`STAB1`](@ref) `ext`."""
set_ext!(value::STAB1, val) = value.ext = val
"""Set [`STAB1`](@ref) `states_types`."""
set_states_types!(value::STAB1, val) = value.states_types = val


text:
module TestAnalysis

using EventTracker
using Plots
using Test

function test_smoke()
    @test EventTracker.summary_dataframe() isa Any
    stks = EventTracker.stacks()
    @test show(devnull, "image/png", plot(stks)) isa Any
end

end  # module


text:
include(joinpath(@__DIR__, "../../src/FLS.jl"))

"""
Semigroups.
"""

p1p2, p2p3, idk, kid = map(Symbol, ["π₁×π₂","π₂×π₃","id×k","k×id"])
semig = FLSinit(@acset LabeledGraph begin
    V = 3
    E = 10
    vlabel = [:s, :s2, :s3]
    elabel = [:k, :π₁, :π₂, :Π₁, :Π₂, :Π₃, p1p2, p2p3, idk, kid]
    src    = [2,  2,   2,   3,    3,   3,   3,   3,    3,   3]
    tgt    = [1,  1,   1,   1,    1,   1,   2,   2,    2,   2]
end );

# s2 is pair
pairconed = @acset LabeledGraph begin
    V = 2
    vlabel = [:s, :s]
end

paircone = Cone(pairconed, :s3, [1=>:π₁, 2=>:π₂])

# s3 is triple
tripconed = @acset LabeledGraph begin
    V = 3
    vlabel = [:s, :s, :s]
end

tripcone = Cone(tripconed, :s3, [1=>:Π₁, 2=>:Π₂, 3=>:Π₃])

semicones = [paircone, tripcone]
semieqs = [
    (:p1p2_p1, [p1p2, :π₁], [:Π₁]),
    (:p1p2_p2, [p1p2, :π₂], [:Π₂]),
    (:p2p3_p1, [p2p3, :π₁], [:Π₂]),
    (:p2p3_p2, [p2p3, :π₂], [:Π₃]),

    (:kid_p1, [kid, :π₁], [p1p2,:k]),
    (:kid_p2, [kid, :π₂], [:Π₃]),

    (:idk_p2, [idk, :π₂], [p2p3,:k]),
    (:idk_p1, [idk, :π₁], [:Π₁]),

    (:assoc, [idk, :k], [kid,:k]),
]
semig = FLS(:semig, semigschema, semigcones, semigeqs);


text:
const outputname = joinpath(@__DIR__, "country-data.jl")

# First, check if currency-data.jl already exists
isfile(outputname) && exit()

using JSON3

const src = "https://pkgstore.datahub.io/core/country-codes/country-codes_json/data/471a2e653140ecdd7243cdcacfd66608/country-codes_json.json"

inputname = joinpath(@__DIR__, "country-codes.json")

const country_list = Dict{String,Tuple}()

function genfile(io)
    for c in country_data
        (cntry = c["ISO3166-1-Alpha-2"]) === nothing && continue
        (nm = c["official_name_en"]) === nothing && continue
        (cd = c["ISO3166-1-numeric"]) === nothing && continue
        ccy_str = c["ISO4217-currency_alphabetic_code"]
        ccy = ccy_str === nothing ? Symbol[] : Symbol.(split(ccy_str,","))
        d = c["Dial"]
        cap = c["Capital"]
        cont = c["Continent"]
        isdev = isequal(c["Developed / Developing Countries"],"Developing")
        reg = c["Region Name"]
        sreg = c["Sub-region Name"]

        country_list[cntry] = (nm,cd,ccy,d,cap,cont,isdev,reg,sreg)
    end
    println(io, "const _country_data = Dict(")
    for (c, val) in country_list
        # println(io, "    :$curr => (Currency{:$curr}, $(val[1]), $(lpad(val[2], 4)), \"$(val[3])\"),")

        print(io, "    :$c => (Country{:$c},",
            "\"$(val[1])\",",
            val[2],",",
            val[3],",",
            "\"$(val[4])\",",
            "\"$(val[5])\",",
            "\"$(val[6])\",",
            val[7],",",
            "\"$(val[8])\",",
            "\"$(val[9])\"),\n")
    end
    println(io, ")\n")

    #     print(io, "    :$(country) => (Country{:$(country)},",
    #         "\"$(official_name)\",",
    #         code,",",
    #         currencies,",",
    #         "\"$(dial)\",",
    #         "\"$(capital)\",",
    #         "\"$(continent)\",",
    #         isdeveloping,",",
    #         "\"$(region)\",",
    #         "\"$(subregion)\"),\n")
    # end
    # println(io, ")\n")
end

const country_data = open(io -> JSON3.read(io), inputname)

open(io -> genfile(io), outputname, "w")


text:
function check_dist_numerical(dist, chn; mean_tol = 0.1, var_atol = 1.0, var_tol = 0.5)
    @testset "numerical" begin
        # Extract values.
        chn_xs = chn[1:2:end, :x, :].value

        # Check means.
        dist_mean = mean(dist)
        mean_shape = size(dist_mean)
        if !all(isnan, dist_mean) && !all(isinf, dist_mean)
            chn_mean = Array(mean(chn_xs, dims=1))
            chn_mean = length(chn_mean) == 1 ?
                chn_mean[1] :
                reshape(chn_mean, mean_shape)
            atol_m = length(chn_mean) > 1 ?
                mean_tol * length(chn_mean) :
                max(mean_tol, mean_tol * chn_mean)
            @test chn_mean ≈ dist_mean atol=atol_m
        end

        # Check variances.
        # var() for Distributions.MatrixDistribution is not defined
        if !(dist isa Distributions.MatrixDistribution)
            # Variance
            dist_var = var(dist)
            var_shape = size(dist_var)
            if !all(isnan, dist_var) && !all(isinf, dist_var)
                chn_var = Array(var(chn_xs, dims=1))
                chn_var = length(chn_var) == 1 ?
                    chn_var[1] :
                    reshape(chn_var, var_shape)
                atol_v = length(chn_mean) > 1 ?
                    mean_tol * length(chn_mean) :
                    max(mean_tol, mean_tol * chn_mean)
                @test chn_mean ≈ dist_mean atol=atol_v
            end
        end
    end
end

# Helper function for numerical tests
function check_numerical(chain,
                        symbols::Vector,
                        exact_vals::Vector;
                        atol=0.2,
                        rtol=0.0)
    for (sym, val) in zip(symbols, exact_vals)
        E = val isa Real ?
            mean(chain[sym].value) :
            vec(mean(chain[sym].value, dims=[1]))
        @info (symbol=sym, exact=val, evaluated=E)
        @test E ≈ val atol=atol rtol=rtol
    end
end

# Wrapper function to quickly check gdemo accuracy.
function check_gdemo(chain; atol=0.2, rtol=0.0)
    check_numerical(chain, [:s, :m], [49/24, 7/6], atol=atol, rtol=rtol)
end

# Wrapper function to check MoGtest.
function check_MoGtest_default(chain; atol=0.2, rtol=0.0)
    check_numerical(chain,
        [:z1, :z2, :z3, :z4, :mu1, :mu2],
        [1.0, 1.0, 2.0, 2.0, 1.0, 4.0],
        atol=atol, rtol=rtol)
end


text:
using TestFailure
using Test

if haskey(ENV, "TEST_SIGNAL")
    run(`kill -s $(ENV["TEST_SIGNAL"]) $(getpid())`)
elseif haskey(ENV, "TEST_EXITCODE")
    exit(parse(Int, ENV["TEST_EXITCODE"]))
end

@testset "TestFailure" begin
    @test false
end


text:
using CUDA

###########################################################################
## helper functions (originally from CuArrays, but since removed there)  ##
###########################################################################
function cudims(n::Integer)
  threads = min(n, 256)
  cld(n,threads), threads
end

cudims(a::AbstractArray) = cudims(length(a))
cudims(a::SymArray) = cudims(a.data)

# COV_EXCL_START
@inline ind2sub_(a::AbstractArray{T,0}, i) where T = ()
@inline ind2sub_(a, i) = Tuple(CartesianIndices(a)[i])

macro cuindex(A)
  quote
    A = $(esc(A))
    i = (blockIdx().x-1) * blockDim().x + threadIdx().x
    i > length(A) && return
    ind2sub_(A, i)
  end
end
# COV_EXCL_STOP
###########################################################################

mygemv!(tA,alpha,A::CuArray,args...) = CUDA.CUBLAS.gemv!(tA,alpha,A,args...)
_contract_middle!(res::CuArray,A,B) = (@tensor res[i,k] = B[i,j,k] * A[j])

# for SymArrays on the GPU, collect should convert to a SymArray on the CPU
# to convert to a "normal" Array, you then have to apply collect again
Base.collect(S::SymArray{Nsyms,T,N,M,datType}) where {Nsyms,T,N,M,datType<:CuArray} = SymArray{Nsyms}(collect(S.data),S.size...)

@generated function cuda_contraction_kernel(res, A, S, SI::SymIndexIter{Nsymres}) where {Nsymres}
    # calculate res[iA1,iA3,iS1,iSm2,iS3] = ∑_iA2 A[iA1,iA2,iA3] * S[iS1,iS2,iS3]
    # where iSm2 = (i1,i2,...,iNsymres) and iS2 = sorted(iA2,i1,i2,i3...,iNsymres)
    code = quote
        I = @cuindex(res)
        iA1,iA3,iS1,iSm2,iS3 = I
        ISm2 = ind2sub_symgrp(SI, iSm2)
        res[I...] = zero(eltype(res))
    end
    for n = 0:Nsymres
        iAstart = n==0 ? 1 : :(ISm2[$n]+1)
        iAend = n<Nsymres ? :(ISm2[$(n+1)]) : :(size(A,2))
        iprev = [:( ISm2[$i] ) for i=1:n]
        ipost = [:( ISm2[$i] ) for i=n+1:Nsymres]
        cc = :( for iA2 = $iAstart:$iAend
                    iS2 = symgrp_sortedsub2ind($(iprev...),iA2,$(ipost...))
                    res[I...] += A[iA1,iA2,iA3]*S[iS1,iS2,iS3]
                end)
        push!(code.args,cc)
    end
    push!(code.args,:(return))
    #display(code)
    :( @inbounds $code )
end

function contract_symindex!(res::CuArray{TU,5}, A::CuArray{T,3}, ::Val{sizeA13unit}, S::CuArray{U,3}, ::Val{sizeS13unit}, ::Val{Nsym}) where {T,U,TU,sizeA13unit,sizeS13unit,Nsym}
    blk, thr = cudims(res)
    SI = SymIndexIter(Nsym-1,size(A,2))
    @cuda blocks=blk threads=thr cuda_contraction_kernel(res,A,S,SI)
end


text:
function randstrarray(pool, N)
  K = length(pool)
  PooledArray(PooledArrays.RefArray(rand(UInt8(1):UInt8(K), N)), pool)
end

function randstrarray1(pool, N)
    K = length(pool)
    PooledArray(PooledArrays.RefArray(rand(1:K, N)), pool)
end

function randstrarrayc(pool, N)
  K = length(pool)
  CategoricalArray(rand(1:K,N) .|> UInt32,CategoricalPool(pool))
end

function randstrarray1c(pool, N)
    K = length(pool)
    CategoricalArray(rand(1:K,N) .|> UInt32,CategoricalPool(pool))
end

function createIndexedTable(N::Int,K::Int)
  pool = [@sprintf "id%03d" k for k in 1:K]
  pool1 = [@sprintf "id%010d" k for k in 1:(N/K)]

  DT = IndexedTable(
    IndexedTables.Columns(
      row_id = [1:N;]
      ),
    IndexedTables.Columns(
      id1 = randstrarray(pool, N),
      id2 = randstrarray(pool, N),
      id3 = randstrarray1(pool1, N),
      id4 = rand(1:K, N),                          # large groups (int)
      id5 = rand(1:K, N),                          # large groups (int)
      id6 = rand(1:Int64(N/K), N),                        # small groups (int)
      v1 =  rand(1:5, N),                          # int in range [1,5]
      v2 =  rand(1:5, N),                          # int in range [1,5]
      v3 =  rand(round.(rand(100)*100,4), N) # numeric e.g. 23.5749
      ));
  return DT
end

function run_juliadb_bench(N::Int = Int64(2e9/8), K::Int = 100)
  # warm up
  createIndexedTable(1,1)

  # create true dataset
  #julia_timings = Dict{Symbol, Float64}()
  timegen = @elapsed DT = createIndexedTable(N, K)
  julia_timings = Dict(
    :gen_syn => timegen,
    :sum1_1 => (@elapsed sumby(DT, :id1, :v1)),
    :sum1_2 => (@elapsed sumby(DT, :id1, :v1)),
    :sum2_1 => (@elapsed groupreduce(x->(x[1], x[2]), x -> x[3], (x, y) -> x+y, zip(column(DT,:id1),column(DT,:id2),column(DT,:v1)))),
    :sum2_2 => (@elapsed groupreduce(x->(x[1], x[2]), x -> x[3], (x, y) -> x+y, zip(column(DT,:id1),column(DT,:id2),column(DT,:v1)))),
    :sum_mean1 => (@elapsed [sumby(DT, :id3, :v1), meanby(DT,:id3, :v3)]),
    :sum_mean2 => (@elapsed [sumby(DT, :id3, :v1), meanby(DT,:id3, :v3)]),
    :mean7_9_by_id4_1 => (@elapsed [meanby(DT, :id4, :v1), meanby(DT, :id4, :v2), meanby(DT, :id4, :v3)]),
    :mean7_9_by_id4_2 => (@elapsed [meanby(DT, :id4, :v1), meanby(DT, :id4, :v2), meanby(DT, :id4, :v3)]),
    :sum7_9_by_id6_1 => (@elapsed groupreduce(
        x->x[1],
        x->(x[2],x[3],x[4]),
        (x,y)->(x[1]+y[1],x[2]+y[2],x[3]+y[3]),
        zip(column(DT,:id6),column(DT,:v1),column(DT,:v2),column(DT,:v3)))),
    :sum7_9_by_id6_2 => (@elapsed groupreduce(
        x->x[1],
        x->(x[2],x[3],x[4]),
        (x,y)->(x[1]+y[1],x[2]+y[2],x[3]+y[3]),
        zip(column(DT,:id6),column(DT,:v1),column(DT,:v2),column(DT,:v3))))
  )

  return julia_timings
end

function run_juliadb_bench_pmap(N::Int = Int64(2e9/8), K::Int = 100)
  # warm up
  createIndexedTable(1,1)

  # create true dataset
  #julia_timings = Dict{Symbol, Float64}()
  timegen = @elapsed DT = createIndexedTable(N, K)
  julia_timings = Dict(
    :gen_syn => timegen,
    :sum1_1 => (@elapsed psumby(DT, :id1, :v1)),
    :sum1_2 => (@elapsed psumby(DT, :id1, :v1)),
    :sum2_1 => (@elapsed groupreduce(x->(x[1], x[2]), x -> x[3], (x, y) -> x+y, zip(column(DT,:id1),column(DT,:id2),column(DT,:v1)))),
    :sum2_2 => (@elapsed groupreduce(x->(x[1], x[2]), x -> x[3], (x, y) -> x+y, zip(column(DT,:id1),column(DT,:id2),column(DT,:v1)))),
    :sum_mean1 => (@elapsed [psumby(DT, :id3, :v1), pmeanby(DT,:id3, :v3)]),
    :sum_mean2 => (@elapsed [psumby(DT, :id3, :v1), pmeanby(DT,:id3, :v3)]),
    :mean7_9_by_id4_1 => (@elapsed [pmeanby(DT, :id4, :v1), pmeanby(DT, :id4, :v2), pmeanby(DT, :id4, :v3)]),
    :mean7_9_by_id4_2 => (@elapsed [pmeanby(DT, :id4, :v1), pmeanby(DT, :id4, :v2), pmeanby(DT, :id4, :v3)]),
    :sum7_9_by_id6_1 => (@elapsed groupreduce(
        x->x[1],
        x->(x[2],x[3],x[4]),
        (x,y)->(x[1]+y[1],x[2]+y[2],x[3]+y[3]),
        zip(column(DT,:id6),column(DT,:v1),column(DT,:v2),column(DT,:v3)))),
    :sum7_9_by_id6_2 => (@elapsed groupreduce(
        x->x[1],
        x->(x[2],x[3],x[4]),
        (x,y)->(x[1]+y[1],x[2]+y[2],x[3]+y[3]),
        zip(column(DT,:id6),column(DT,:v1),column(DT,:v2),column(DT,:v3))))
  )

  return julia_timings
end

# function run_sac_bench(N::Int = Int64(2e9/8), K::Int = 100)
#   # create true dataset
#   #julia_timings = Dict{Symbol, Float64}()
#   timegen = @elapsed DT = createIndexedTable(N, K)
#   sac_timings = Dict(
#     :gen_syn => timegen,
#     :sum1_1 => (@elapsed groupreduce(x->x[1], x->x[2], (x,y) -> x+y, zip(column(DT,:id1),column(DT,:v1)))),
#     :sum1_2 => (@elapsed groupreduce(x->x[1], x->x[2], (x,y) -> x+y, zip(column(DT,:id1),column(DT,:v1)))),
#     :sum2_1 => (@elapsed groupreduce(x->(x[1], x[2]), x -> x[3], (x, y) -> x+y, zip(column(DT,:id1),column(DT,:id2),column(DT,:v1)))),
#     :sum2_2 => (@elapsed groupreduce(x->(x[1], x[2]), x -> x[3], (x, y) -> x+y, zip(column(DT,:id1),column(DT,:id2),column(DT,:v1)))),
#     :sum_mean1 => (@elapsed groupreduce(x->x[1], x -> (x[2],x[3],1), (x, y) -> (x[1]+y[1],x[2]+y[2],x[3]+y3]), zip(column(DT,:id3),column(DT,:v1),column(DT,:v3)))),
#     :sum_mean2 => (@elapsed groupreduce(x->x[1], x -> (x[2],x[3],1), (x, y) -> (x[1]+y[1],x[2]+y[2],x[3]+y3]), zip(column(DT,:id3),column(DT,:v1),column(DT,:v3)))),
#     :mean7_9_by_id4_1 => (@elapsed [meanby(DT, :id4, :v1), meanby(DT, :id4, :v2), meanby(DT, :id4, :v3)]),
#     :mean7_9_by_id4_2 => (@elapsed [meanby(DT, :id4, :v1), meanby(DT, :id4, :v2), meanby(DT, :id4, :v3)]),
#     :sum7_9_by_id6_1 => (@elapsed [sumby(DT, :id6, :v1), sumby(DT, :id6, :v2), sumby(DT, :id6, :v3)]),
#     :sum7_9_by_id6_2 => (@elapsed [sumby(DT, :id6, :v1), sumby(DT, :id6, :v2), sumby(DT, :id6, :v3)])
#   )
#
#   return sac_timings
# end


text:
function Jdouble(x::Vector{Float64},beta::Vector{Float64},rho::Float64)
    J = zeros(6,6);

    J[1,1] = 1.;
    J[2,1] = rho.*x[1]./(x[2].^2);
    J[3,1] = 1.;
    J[4,1] = 1./x[2];

    J[2,2] = (0.5.*beta[1].*x[2].^-0.5 .- rho.*x[1].^2.*x[2].^-3);
    J[3,2] = J[2,2];
    J[4,2] = -x[1].*x[2].^-2 .+ sqrt.(0.5.*beta[1]./rho).*x[2].^-0.75;

    J[1,3] = -1.;
    J[2,3] = -rho.*x[3]./(x[4].^2);
    J[5,3] = 1./x[4];

    J[2,4] = (-0.5.*beta[2].*x[4].^-0.5 .+ rho.*x[3].^2.*x[4].^-3);
    J[5,4] = -x[3].*x[4].^-2 .- sqrt.(0.5.*beta[2]./rho).*x[4].^-0.75;

    J[1,5] = -1.;
    J[3,5] = -rho.*x[5]./(x[6].^2);
    J[6,5] = 1./x[6];

    J[3,6] = (-0.5.*beta[3].*x[6].^-0.5 .+ rho.*x[5].^2.*x[6].^-3);
    J[6,6] = -x[5].*x[6].^-2 .- sqrt.(0.5.*beta[3]./rho).*x[6].^-0.75;

    return J
end

max_stars_repo_name: dcanuto/CardioModeling

max_issues_repo_name: dcanuto/CardioModeling

max_forks_repo_name: dcanuto/CardioModeling


text:
using Singular
using AbstractAlgebra

"""
    @default ex, fs

For each `f` in `fs` declare a method `f(args...; kwargs...) = f(ex, args...; kwargs)`.
"""

macro default(ex, fs)
    ex = esc(ex)
    fs = isexpr(fs, :tuple) ? map(esc, fs.args) : [esc(fs)]
    quote
        $([:($f(args...; kwargs...) = $f($ex, args...; kwargs...)) for f in fs]...)
    end
end

#.. Helpers for BitVectors

function bitarray_to_int(arr::BitArray,
                         B::Type{U}) where {U <: Integer}
    arr = reverse(arr)
    sum(((i, x),) -> B(x) << ((i-1) * sizeof(x)), enumerate(arr.chunks))
end

function even_partition(i, nums)
    fl = typeof(i)(floor(i / nums))
    part = repeat([fl], nums)
    while sum(part) != i
        ind = rand(1:nums)
        part[ind] += 1
    end
    part
end

function even_between(a, b, nums)
    part = even_partition(b - a, nums)
    part_2 = similar(part)
    curr = a
    for (i, p) in enumerate(part)
        part_2[i] = curr
        curr += p
    end
    part_2
end

function mac_bound(I::Vector{P}) where {P <: AbstractAlgebra.MPolyElem}
    I_sorted = sort(I, by = p -> Singular.total_degree(p), rev = true)
    l = min(Singular.nvars(parent(first(I))) + 1, length(I))
    sum([total_degree(I_sorted[j]) for j in 1:l]) - l + 1
end


text:
"""
    IndicatorSOS1Bridge{T, BC <: MOI.AbstractScalarSet}

The `IndicatorSOS1Bridge` replaces an indicator constraint of the following form:
``z \\in \\mathbb{B}, z == 1 \\implies f(x) \\leq b`` with a SOS1 constraint:
``z \\in \\mathbb{B}, w \\leq 0, f(x) + w \\leq b, SOS1(w, z)``.
`GreaterThan` constraints are handled in a symmetric way:
``z \\in \\mathbb{B}, z == 1 \\implies f(x) \\geq b`` is reformulated as:
``z \\in \\mathbb{B}, w \\geq 0, f(x) + w \\geq b, SOS1(w, z)``.
Other scalar sets are handled without a bound constraint:
``z \\in \\mathbb{B}, z == 1 \\implies f(x) == b`` is reformulated as:
``z \\in \\mathbb{B}, w \\text{ free}, f(x) + w == b, SOS1(w, z)``.

If `BC !<: Union{LessThan, GreaterThan}`, `bound_constraint_index` is `nothing`.
"""
struct IndicatorSOS1Bridge{
    T,
    BC<:MOI.AbstractScalarSet,
    MaybeBC<:Union{MOI.ConstraintIndex{MOI.SingleVariable,BC},Nothing},
} <: AbstractBridge
    w_variable_index::MOI.VariableIndex
    z_variable_index::MOI.VariableIndex
    affine_func::MOI.ScalarAffineFunction{T}
    bound_constraint_index::MaybeBC
    sos_constraint_index::MOI.ConstraintIndex{MOI.VectorOfVariables,MOI.SOS1{T}}
    linear_constraint_index::MOI.ConstraintIndex{MOI.ScalarAffineFunction{T},BC}
end

function bridge_constraint(
    ::Type{IndicatorSOS1Bridge{T,BC,MaybeBC}},
    model::MOI.ModelLike,
    f::MOI.VectorAffineFunction{T},
    s::MOI.IndicatorSet{MOI.ACTIVATE_ON_ONE,BC},
) where {T<:Real,BC,MaybeBC}
    f_scalars = MOIU.eachscalar(f)
    (w, bound_constraint) = _add_bound_constraint!(model, BC)
    z = convert(MOI.SingleVariable, f_scalars[1]).variable
    sos_vector = MOI.VectorOfVariables([w, z])
    sos_constraint =
        MOI.add_constraint(model, sos_vector, MOI.SOS1{T}([0.4, 0.6]))
    affine_func = f_scalars[2]
    affine_expr = MOIU.operate(+, T, affine_func, MOI.SingleVariable(w))
    linear_constraint = MOI.add_constraint(model, affine_expr, s.set)
    return IndicatorSOS1Bridge{T,BC,MaybeBC}(
        w,
        z,
        affine_func,
        bound_constraint,
        sos_constraint,
        linear_constraint,
    )
end

function _add_bound_constraint!(
    model::MOI.ModelLike,
    ::Type{BC},
) where {T<:Real,BC<:Union{MOI.LessThan{T},MOI.GreaterThan{T}}}
    return MOI.add_constrained_variable(model, BC(zero(T)))
end

function _add_bound_constraint!(
    model::MOI.ModelLike,
    ::Type{<:MOI.AbstractScalarSet},
)
    return (MOI.add_variable(model), nothing)
end

function MOI.supports_constraint(
    ::Type{<:IndicatorSOS1Bridge},
    ::Type{<:MOI.AbstractVectorFunction},
    ::Type{<:MOI.IndicatorSet{MOI.ACTIVATE_ON_ONE,<:MOI.AbstractScalarSet}},
)
    return true
end

function MOI.get(
    model::MOI.ModelLike,
    attr::MOI.ConstraintSet,
    b::IndicatorSOS1Bridge,
)
    return MOI.IndicatorSet{MOI.ACTIVATE_ON_ONE}(
        MOI.get(model, attr, b.linear_constraint_index),
    )
end

function MOI.get(
    model::MOI.ModelLike,
    attr::MOI.ConstraintFunction,
    b::IndicatorSOS1Bridge{T},
) where {T}
    z = b.z_variable_index
    terms = [MOI.VectorAffineTerm(1, MOI.ScalarAffineTerm(one(T), z))]
    for affine_term in b.affine_func.terms
        push!(terms, MOI.VectorAffineTerm(2, affine_term))
    end
    return MOI.VectorAffineFunction(terms, [zero(T), b.affine_func.constant])
end

function MOI.delete(model::MOI.ModelLike, bridge::IndicatorSOS1Bridge)
    if bridge.bound_constraint_index !== nothing
        MOI.delete(model, bridge.bound_constraint_index)
    end
    MOI.delete(model, bridge.sos_constraint_index)
    MOI.delete(model, bridge.linear_constraint_index)
    MOI.delete(model, bridge.w_variable_index)
    return
end

function MOIB.added_constrained_variable_types(
    ::Type{<:IndicatorSOS1Bridge{T,BC}},
) where {T,BC<:Union{MOI.LessThan{T},MOI.GreaterThan{T}}}
    return [(BC,)]
end

function MOIB.added_constrained_variable_types(
    ::Type{<:IndicatorSOS1Bridge{T,BC}},
) where {T,BC}
    return []
end

function MOIB.added_constraint_types(
    ::Type{<:IndicatorSOS1Bridge{T,BC}},
) where {T,BC<:Union{MOI.LessThan{T},MOI.GreaterThan{T}}}
    return [
        (MOI.VectorOfVariables, MOI.SOS1{T}),
        (MOI.ScalarAffineFunction{T}, BC),
    ]
end

function MOIB.added_constraint_types(
    ::Type{<:IndicatorSOS1Bridge{T,S}},
) where {T,S<:MOI.AbstractScalarSet}
    return [
        (MOI.VectorOfVariables, MOI.SOS1{T}),
        (MOI.ScalarAffineFunction{T}, S),
    ]
end

function concrete_bridge_type(
    ::Type{<:IndicatorSOS1Bridge{T}},
    ::Type{<:MOI.AbstractVectorFunction},
    ::Type{MOI.IndicatorSet{MOI.ACTIVATE_ON_ONE,S}},
) where {T,S<:Union{MOI.LessThan,MOI.GreaterThan}}
    return IndicatorSOS1Bridge{T,S,MOI.ConstraintIndex{MOI.SingleVariable,S}}
end

function concrete_bridge_type(
    ::Type{<:IndicatorSOS1Bridge{T}},
    ::Type{<:MOI.AbstractVectorFunction},
    ::Type{MOI.IndicatorSet{MOI.ACTIVATE_ON_ONE,S}},
) where {T,S<:MOI.AbstractScalarSet}
    return IndicatorSOS1Bridge{T,S,Nothing}
end

# Attributes, Bridge acting as a model

function MOI.get(::IndicatorSOS1Bridge, ::MOI.NumberOfVariables)
    return 1
end

function MOI.get(b::IndicatorSOS1Bridge, ::MOI.ListOfVariableIndices)
    return [b.w_variable_index]
end

function MOI.get(
    ::IndicatorSOS1Bridge{T,BC,Nothing},
    ::MOI.NumberOfConstraints{MOI.SingleVariable,BC},
) where {T,BC}
    return 0
end

function MOI.get(
    ::IndicatorSOS1Bridge{T,BC,CI},
    ::MOI.NumberOfConstraints{MOI.SingleVariable,BC},
) where {T,BC,CI<:MOI.ConstraintIndex{MOI.SingleVariable,BC}}
    return 1
end

function MOI.get(
    ::IndicatorSOS1Bridge,
    ::MOI.NumberOfConstraints{MOI.VectorOfVariables,<:MOI.SOS1},
)
    return 1
end

function MOI.get(
    ::IndicatorSOS1Bridge{T,BC},
    ::MOI.NumberOfConstraints{MOI.ScalarAffineFunction{T},BC},
) where {T,BC,CI<:MOI.ConstraintIndex{MOI.SingleVariable,BC}}
    return 1
end

function MOI.get(
    b::IndicatorSOS1Bridge{T,BC,CI},
    ::MOI.ListOfConstraintIndices{MOI.SingleVariable,BC},
) where {T,BC,CI<:MOI.ConstraintIndex}
    return [b.bound_constraint_index]
end

function MOI.get(
    ::IndicatorSOS1Bridge{T,BC,Nothing},
    ::MOI.ListOfConstraintIndices{MOI.SingleVariable,BC},
) where {T,BC}
    return MOI.ConstraintIndex{MOI.SingleVariable,BC}[]
end

function MOI.get(
    b::IndicatorSOS1Bridge{T},
    ::MOI.ListOfConstraintIndices{MOI.VectorOfVariables,<:MOI.SOS1},
) where {T}
    return [b.sos_constraint_index]
end

function MOI.get(
    b::IndicatorSOS1Bridge{T,BC},
    ::MOI.ListOfConstraintIndices{MOI.ScalarAffineFunction{T},BC},
) where {T,BC}
    return [b.linear_constraint_index]
end

function MOI.supports(
    ::MOI.ModelLike,
    ::MOI.ConstraintPrimalStart,
    ::Type{<:IndicatorSOS1Bridge},
)
    return true
end

function MOI.get(
    model::MOI.ModelLike,
    attr::MOI.ConstraintPrimal,
    bridge::IndicatorSOS1Bridge,
)
    zvalue = MOI.get(model, MOI.VariablePrimal(attr.N), bridge.z_variable_index)
    wvalue = MOI.get(model, MOI.VariablePrimal(attr.N), bridge.w_variable_index)
    lin_primal_start = MOI.get(model, attr, bridge.linear_constraint_index)
    return [zvalue, lin_primal_start - wvalue]
end

function MOI.get(
    model::MOI.ModelLike,
    attr::MOI.ConstraintPrimalStart,
    bridge::IndicatorSOS1Bridge,
)
    zstart = MOI.get(model, MOI.VariablePrimalStart(), bridge.z_variable_index)
    wstart = MOI.get(model, MOI.VariablePrimalStart(), bridge.w_variable_index)
    lin_primal_start = MOI.get(model, attr, bridge.linear_constraint_index)
    return [zstart, lin_primal_start - wstart]
end

function MOI.set(
    model::MOI.ModelLike,
    attr::MOI.ConstraintPrimalStart,
    bridge::IndicatorSOS1Bridge{T},
    value,
) where {T}
    zvalue = value[1]
    lin_start = value[2]
    MOI.set(model, MOI.VariablePrimalStart(), bridge.z_variable_index, zvalue)
    wstart = MOI.get(model, MOI.VariablePrimalStart(), bridge.w_variable_index)
    wstart = wstart === nothing ? zero(T) : wstart
    return MOI.set(
        model,
        attr,
        bridge.linear_constraint_index,
        lin_start + wstart,
    )
end


text:
import Base.==, Base.+

export ==, EqConstraint, +

AffineOrConst = Union{AffineExpr, Constant}

type EqConstraint
  head::Symbol
  lhs::AffineOrConst
  rhs::AffineOrConst
  canon_form::AffineExpr

  function EqConstraint(lhs::AffineOrConst, rhs::AffineOrConst)
    if lhs.head == :Constant && rhs.head == :Constant
      error("Equality constraints between two constants are not allowed")
    end
    if lhs.size != rhs.size && lhs.size != (1, 1) && rhs.size != (1, 1)
      error("LHS and RHS of constraints muust have the same size, or one needs to be scalar")
    end
    this = new(:(==), lhs, rhs, lhs - rhs)
    return this
  end
end

==(lhs::AffineOrConst, rhs::AffineOrConst) = EqConstraint(lhs, rhs)
==(lhs::Numeric, rhs::AffineExpr) = ==(convert(Constant, lhs), rhs)
==(lhs::AffineExpr, rhs::Numeric) = ==(lhs, convert(Constant, rhs))

+(constraints::Array{EqConstraint, 1}, new_constraints::Array{EqConstraint, 1}) = append!(constraints, new_constraints)
+(constraints::Array{EqConstraint, 1}, new_constraint::EqConstraint) = push!(constraints, new_constraint)
+(constraints::Array{Any, 1}, new_constraints::Array{EqConstraint, 1}) =
  if length(constraints) == 0
    new_constraints
  else
    error("Cannot append equality constraint to arbitrary array.")
  end
+(constraints::Array{Any, 1}, new_constraint::EqConstraint) =
  if length(constraints) == 0
    [new_constraint]
  else
    error("Cannot append equality constraint to arbitrary array.")
  end