text:
#################### Core Model Functionality ####################

#################### Constructors ####################

"""
    Model(; iter::Integer=0, burnin::Integer=0,
      samplers::Vector{Sampler}=Sampler[], nodes...)

Construct a `Model` object that defines a model for MCMC simulation.

Returns a `Model` type object.

* `iter`: current iteration of the MCMC simulation.

* `burnin`: number of initial draws to be discarded as a burn-in sequence to allow for convergence.

* `samplers`: block-specific sampling functions.

* `nodes...`: arbitrary number of user-specified arguments defining logical and stochastic nodes in the model. Argument values must be `Logical` or `Stochastic` type objects. Their names in the model will be taken from the argument names.

"""
function Model(; iter::Integer=0, burnin::Integer=0,
               samplers::Vector{Sampler}=Sampler[], nodes...)

  nodedict = Dict{Symbol, Any}()
  for (key, value) in nodes
    isa(value, AbstractDependent) || throw(ArgumentError("nodes are not all Dependent types"))
    node = deepcopy(value)
    isa(key, Symbol) || throw(ArgumentError("Something is wrong here"))
    node.symbol = key
    nodedict[key] = node
  end
  m = Model(nodedict, Sampler[], ModelState[], iter, burnin, false, false, -Inf64)
  dag = ModelGraph(m)
  dependentkeys = keys(m, :dependent)
  terminalkeys = keys(m, :stochastic)
  for v in vertices(dag.graph)
    vkey = dag.keys[v]
    if vkey in dependentkeys
      m[vkey].targets = intersect(dependentkeys,
                                  gettargets(dag, v, terminalkeys))
    end
  end
  setsamplers!(m, samplers)
end

#################### Indexing ####################

Base.getindex(m::Model, nodekey::Symbol) = m.nodes[nodekey]

function Base.setindex!(m::Model, value, nodekey::Symbol)
  node = m[nodekey]
  if isa(node, AbstractDependent)
    node.value = value
  else
    m.nodes[nodekey] = convert(typeof(node), value)
  end
end

function Base.setindex!(m::Model, values::Dict, nodekeys::Vector{Symbol})
  for key in nodekeys
    m[key] = values[key]
  end
end

function Base.setindex!(m::Model, value, nodekeys::Vector{Symbol})
  length(nodekeys) == 1 || throw(BoundsError())
  m[first(nodekeys)] = value
end

Base.keys(m::Model) = collect(keys(m.nodes))
"""
    Base.keys(m::Model, ntype::Symbol, at...)

Extract the symbols (keys) for all existing nodes or for nodes of a specified type.

* `m` : model containing the nodes of interest.

* `ntype` : type of nodes to return. Options are
  * `:all` : all input, logical, and stochastic model nodes.

  * `:assigned` : nodes that have been assigned values.

  * `:block` : stochastic nodes being updated by the sampling block(s) `at::Integer=0` (default: all blocks).

  * `:dependent` : logical and stochastic (dependent) nodes in topologically sorted order.

  * `:independent` or `:input` : input (independent) nodes.

  * `:logical` : logical nodes.

  * `:monitor` : stochastic nodes being monitored in MCMC sampler output.

  * `:output` : stochastic nodes upon which no other stochastic nodes depend.

  * `:source` : nodes upon which the node `at::Symbol` or vector of nodes `at::Vector{Symbol}` depends.

  * `:stochastic` : stochastic nodes.

  * `:target` : topologically sorted nodes that depend on the sampling block(s) `at::Integer=0` (default: all blocks), node `at::Symbol` , or vector of nodes `at::Vector{Symbol}` .

* `at...` : additional positional arguments to be passed to the `ntype` options, as described above.
"""
function Base.keys(m::Model, ntype::Symbol, at...)
  ntype == :block       ? keys_block(m, at...) :
  ntype == :all         ? keys_all(m) :
  ntype == :assigned    ? keys_assigned(m) :
  ntype == :dependent   ? keys_dependent(m) :
  ntype == :independent ? keys_independent(m) :
  ntype == :input       ? keys_independent(m) :
  ntype == :logical     ? keys_logical(m) :
  ntype == :monitor     ? keys_monitor(m) :
  ntype == :output      ? keys_output(m) :
  ntype == :source      ? keys_source(m, at...) :
  ntype == :stochastic  ? keys_stochastic(m) :
  ntype == :target      ? keys_target(m, at...) :
    throw(ArgumentError("unsupported node type $ntype"))
end

function keys_all(m::Model)::Array{Symbol}
  values = Symbol[]
  for key in keys(m)
    node = m[key]
    if isa(node, AbstractDependent)
      push!(values, key)
      append!(values, node.sources)
    end
  end
  unique(values)
end

function keys_assigned(m::Model)::Array{Symbol}
  if m.hasinits
    values = keys(m)
  else
    values = Symbol[]
    for key in keys(m)
      if !isa(m[key], AbstractDependent)
        push!(values, key)
      end
    end
  end
  values
end

function keys_block(m::Model, block::Integer=0)::Array{Symbol}
  block == 0 ? keys_block0(m) : m.samplers[block].params
end

function keys_block0(m::Model)::Array{Symbol}
  values = Symbol[]
  for sampler in m.samplers
    append!(values, sampler.params)
  end
  unique(values)
end

function keys_dependent(m::Model)::Array{Symbol}
  values = Symbol[]
  for key in keys(m)
    if isa(m[key], AbstractDependent)
      push!(values, key)
    end
  end
  intersect(tsort(m), values)
end

function keys_independent(m::Model)::Array{Symbol}
  deps = Symbol[]
  for key in keys(m)
    if isa(m[key], AbstractDependent)
      push!(deps, key)
    end
  end
  setdiff(keys(m, :all), deps)
end

function keys_logical(m::Model)::Array{Symbol}
  values = Symbol[]
  for key in keys(m)
    if isa(m[key], AbstractLogical) || isa(m[key], TreeLogical)
      push!(values, key)
    end
  end
  values
end

function keys_monitor(m::Model)::Array{Symbol}
  values = Symbol[]
  for key in keys(m)
    node = m[key]
    if isa(node, AbstractDependent) && !isempty(node.monitor)
      push!(values, key)
    end
  end
  values
end

function keys_output(m::Model)::Array{Symbol}
  values = Symbol[]
  dag = ModelGraph(m)
  for v in vertices(dag.graph)
    vkey = dag.keys[v]
    if isa(m[vkey], AbstractStochastic) && !any_stochastic(dag, v, m)
      push!(values, vkey)
    end
  end
  values
end

keys_source(m::Model, nodekey::Symbol)::Array{Symbol} = m[nodekey].sources

function keys_source(m::Model, nodekeys::Vector{Symbol})::Array{Symbol}
  values = Symbol[]
  for key in nodekeys
    append!(values, m[key].sources)
  end
  unique(values)
end

function keys_stochastic(m::Model)::Array{Symbol}
  values = Symbol[]
  for key in keys(m)
    if isa(m[key], AbstractStochastic) || isa(m[key], TreeStochastic)
      push!(values, key)
    end
  end
  values
end

function keys_target(m::Model, block::Integer=0)::Array{Symbol}
  block == 0 ? keys_target0(m) : m.samplers[block].targets
end

function keys_target0(m::Model)::Array{Symbol}
  values = Symbol[]
  for sampler in m.samplers
    append!(values, sampler.targets)
  end
  intersect(keys(m, :dependent), values)
end

keys_target(m::Model, nodekey::Symbol)::Array{Symbol} = m[nodekey].targets

function keys_target(m::Model, nodekeys::Vector{Symbol})::Array{Symbol}
  values = Symbol[]
  for key in nodekeys
    append!(values, m[key].targets)
  end
  intersect(keys(m, :dependent), values)
end

#################### Display ####################
"""
    Base.show(io::IO, m::Model)

Write a text representation of the model, nodes, and attributes to the current output stream.
"""
function Base.show(io::IO, m::Model)
  showf(io, m, Base.show)
end
"""
    showall(io::IO, m::Model)

Write a verbose text representation of the model, nodes, and attributes to the current output stream.
"""
function showall(io::IO, m::Model)
  showf(io, m, Base.showall)
end

function showf(io::IO, m::Model, f::Function)
  print(io, "Object of type \"$(summary(m))\"\n")
  width = displaysize()[2] - 1
  for node in keys(m)
    print(io, string("-"^width, "\n", node, ":\n"))
    f(io, m[node])
    println(io)
  end
end

#################### Auxiliary Functions ####################

function names(m::Model, monitoronly::Bool)
  values = AbstractString[]
  for key in keys(m, :dependent)
    if monitoronly
      if !isempty(m[key].monitor)
        nodenames = names(m, key)
        append!(values, nodenames)
      end
    else
      nodenames = names(m, key)
      append!(values, nodenames)
    end
    #  end
    #v = monitoronly ? nodenames[m[key].monitor] : nodenames
    #append!(values, v)
  end
  values
end

function names(m::Model, nodekey::Symbol)
  node = m[nodekey]

  unlist(node, names(node))
end

function names(m::Model, nodekeys::Vector{Symbol})
  values = AbstractString[]
  for key in nodekeys
    append!(values, names(m, key))
  end
  values
end

function update_likelihood!(m::Model, likelihood::Float64)
  m.likelihood = likelihood
  m
end


text:
using Plots
using MaxMinFilters
using Random

Random.seed!(0)
len = 300
x = randn(len)
x[1] = 0;
for i = 1+1:len
    x[i] = -(0.5 + x[i-1]*0.8 + x[i]*0.2)
end

w = 5
mx, mn = movmaxmin(x, w)
xrange = movrange(x, w)

plot(x, label = "x")
plot!(mx, label = "mx")
plot!(mn, label = "mn")
plot!(xrange, label = "range")

#png("plot1.png")

envelope = movenvelope(x, w)
plot(x, label = "x")
plot!(envelope[w-1:end], label = "envelope") # with delay fixed

#png("plot2.png")

ext: jl


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

"""
    generate_data_set(data_set_publication_date, data_set_type, destination_s3_bucket_name, role_name_arn, sns_topic_arn)
    generate_data_set(data_set_publication_date, data_set_type, destination_s3_bucket_name, role_name_arn, sns_topic_arn, params::Dict{String,<:Any})

Given a data set type and data set publication date, asynchronously publishes the requested
data set to the specified S3 bucket and notifies the specified SNS topic once the data is
available. Returns a unique request identifier that can be used to correlate requests with
notifications from the SNS topic. Data sets will be published in comma-separated values
(CSV) format with the file name {data_set_type}_YYYY-MM-DD.csv. If a file with the same
name already exists (e.g. if the same data set is requested twice), the original file will
be overwritten by the new file. Requires a Role with an attached permissions policy
providing Allow permissions for the following actions: s3:PutObject, s3:GetBucketLocation,
sns:GetTopicAttributes, sns:Publish, iam:GetRolePolicy.

# Arguments
- `data_set_publication_date`: The date a data set was published. For daily data sets,
  provide a date with day-level granularity for the desired day. For monthly data sets except
  those with prefix disbursed_amount, provide a date with month-level granularity for the
  desired month (the day value will be ignored). For data sets with prefix disbursed_amount,
  provide a date with day-level granularity for the desired day. For these data sets we will
  look backwards in time over the range of 31 days until the first data set is found (the
  latest one).
- `data_set_type`: The desired data set type.
  customer_subscriber_hourly_monthly_subscriptions From 2017-09-15 to present: Available
  daily by 24:00 UTC.   customer_subscriber_annual_subscriptions From 2017-09-15 to present:
  Available daily by 24:00 UTC.   daily_business_usage_by_instance_type From 2017-09-15 to
  present: Available daily by 24:00 UTC.   daily_business_fees From 2017-09-15 to present:
  Available daily by 24:00 UTC.   daily_business_free_trial_conversions From 2017-09-15 to
  present: Available daily by 24:00 UTC.   daily_business_new_instances From 2017-09-15 to
  present: Available daily by 24:00 UTC.   daily_business_new_product_subscribers From
  2017-09-15 to present: Available daily by 24:00 UTC.
  daily_business_canceled_product_subscribers From 2017-09-15 to present: Available daily by
  24:00 UTC.   monthly_revenue_billing_and_revenue_data From 2017-09-15 to present: Available
  monthly on the 15th day of the month by 24:00 UTC. Data includes metered transactions (e.g.
  hourly) from one month prior.   monthly_revenue_annual_subscriptions From 2017-09-15 to
  present: Available monthly on the 15th day of the month by 24:00 UTC. Data includes
  up-front software charges (e.g. annual) from one month prior.
  monthly_revenue_field_demonstration_usage From 2018-03-15 to present: Available monthly on
  the 15th day of the month by 24:00 UTC.   monthly_revenue_flexible_payment_schedule From
  2018-11-15 to present: Available monthly on the 15th day of the month by 24:00 UTC.
  disbursed_amount_by_product From 2017-09-15 to present: Available every 30 days by 24:00
  UTC.   disbursed_amount_by_instance_hours From 2017-09-15 to present: Available every 30
  days by 24:00 UTC.   disbursed_amount_by_customer_geo From 2017-09-15 to present: Available
  every 30 days by 24:00 UTC.   disbursed_amount_by_age_of_uncollected_funds From 2017-09-15
  to present: Available every 30 days by 24:00 UTC.
  disbursed_amount_by_age_of_disbursed_funds From 2017-09-15 to present: Available every 30
  days by 24:00 UTC.   disbursed_amount_by_age_of_past_due_funds From 2018-04-07 to present:
  Available every 30 days by 24:00 UTC.   disbursed_amount_by_uncollected_funds_breakdown
  From 2019-10-04 to present: Available every 30 days by 24:00 UTC.
  sales_compensation_billed_revenue From 2017-09-15 to present: Available monthly on the 15th
  day of the month by 24:00 UTC. Data includes metered transactions (e.g. hourly) from one
  month prior, and up-front software charges (e.g. annual) from one month prior.
  us_sales_and_use_tax_records From 2017-09-15 to present: Available monthly on the 15th day
  of the month by 24:00 UTC.   disbursed_amount_by_product_with_uncollected_funds This data
  set is deprecated. Download related reports from AMMP instead!
  customer_profile_by_industry This data set is deprecated. Download related reports from
  AMMP instead!   customer_profile_by_revenue This data set is deprecated. Download related
  reports from AMMP instead!   customer_profile_by_geography This data set is deprecated.
  Download related reports from AMMP instead!
- `destination_s3_bucket_name`: The name (friendly name, not ARN) of the destination S3
  bucket.
- `role_name_arn`: The Amazon Resource Name (ARN) of the Role with an attached permissions
  policy to interact with the provided AWS services.
- `sns_topic_arn`: Amazon Resource Name (ARN) for the SNS Topic that will be notified when
  the data set has been published or if an error has occurred.

# Optional Parameters
Optional parameters can be passed as a `params::Dict{String,<:Any}`. Valid keys are:
- `"customerDefinedValues"`: (Optional) Key-value pairs which will be returned, unmodified,
  in the Amazon SNS notification message and the data set metadata file. These key-value
  pairs can be used to correlated responses with tracking information from other systems.
- `"destinationS3Prefix"`: (Optional) The desired S3 prefix for the published data set,
  similar to a directory path in standard file systems. For example, if given the bucket name
  \"mybucket\" and the prefix \"myprefix/mydatasets\", the output file \"outputfile\" would
  be published to \"s3://mybucket/myprefix/mydatasets/outputfile\". If the prefix directory
  structure does not exist, it will be created. If no prefix is provided, the data set will
  be published to the S3 bucket root.
"""
function generate_data_set(
    dataSetPublicationDate,
    dataSetType,
    destinationS3BucketName,
    roleNameArn,
    snsTopicArn;
    aws_config::AbstractAWSConfig=global_aws_config(),
)
    return marketplace_commerce_analytics(
        "GenerateDataSet",
        Dict{String,Any}(
            "dataSetPublicationDate" => dataSetPublicationDate,
            "dataSetType" => dataSetType,
            "destinationS3BucketName" => destinationS3BucketName,
            "roleNameArn" => roleNameArn,
            "snsTopicArn" => snsTopicArn,
        );
        aws_config=aws_config,
        feature_set=SERVICE_FEATURE_SET,
    )
end
function generate_data_set(
    dataSetPublicationDate,
    dataSetType,
    destinationS3BucketName,
    roleNameArn,
    snsTopicArn,
    params::AbstractDict{String};
    aws_config::AbstractAWSConfig=global_aws_config(),
)
    return marketplace_commerce_analytics(
        "GenerateDataSet",
        Dict{String,Any}(
            mergewith(
                _merge,
                Dict{String,Any}(
                    "dataSetPublicationDate" => dataSetPublicationDate,
                    "dataSetType" => dataSetType,
                    "destinationS3BucketName" => destinationS3BucketName,
                    "roleNameArn" => roleNameArn,
                    "snsTopicArn" => snsTopicArn,
                ),
                params,
            ),
        );
        aws_config=aws_config,
        feature_set=SERVICE_FEATURE_SET,
    )
end

"""
    start_support_data_export(data_set_type, destination_s3_bucket_name, from_date, role_name_arn, sns_topic_arn)
    start_support_data_export(data_set_type, destination_s3_bucket_name, from_date, role_name_arn, sns_topic_arn, params::Dict{String,<:Any})

Given a data set type and a from date, asynchronously publishes the requested customer
support data to the specified S3 bucket and notifies the specified SNS topic once the data
is available. Returns a unique request identifier that can be used to correlate requests
with notifications from the SNS topic. Data sets will be published in comma-separated
values (CSV) format with the file name {data_set_type}_YYYY-MM-DD'T'HH-mm-ss'Z'.csv. If a
file with the same name already exists (e.g. if the same data set is requested twice), the
original file will be overwritten by the new file. Requires a Role with an attached
permissions policy providing Allow permissions for the following actions: s3:PutObject,
s3:GetBucketLocation, sns:GetTopicAttributes, sns:Publish, iam:GetRolePolicy.

# Arguments
- `data_set_type`:  Specifies the data set type to be written to the output csv file. The
  data set types customer_support_contacts_data and test_customer_support_contacts_data both
  result in a csv file containing the following fields: Product Id, Product Code, Customer
  Guid, Subscription Guid, Subscription Start Date, Organization, AWS Account Id, Given Name,
  Surname, Telephone Number, Email, Title, Country Code, ZIP Code, Operation Type, and
  Operation Time.    customer_support_contacts_data Customer support contact data. The data
  set will contain all changes (Creates, Updates, and Deletes) to customer support contact
  data from the date specified in the from_date parameter.
  test_customer_support_contacts_data An example data set containing static test data in the
  same format as customer_support_contacts_data
- `destination_s3_bucket_name`: The name (friendly name, not ARN) of the destination S3
  bucket.
- `from_date`: The start date from which to retrieve the data set in UTC. This parameter
  only affects the customer_support_contacts_data data set type.
- `role_name_arn`: The Amazon Resource Name (ARN) of the Role with an attached permissions
  policy to interact with the provided AWS services.
- `sns_topic_arn`: Amazon Resource Name (ARN) for the SNS Topic that will be notified when
  the data set has been published or if an error has occurred.

# Optional Parameters
Optional parameters can be passed as a `params::Dict{String,<:Any}`. Valid keys are:
- `"customerDefinedValues"`: (Optional) Key-value pairs which will be returned, unmodified,
  in the Amazon SNS notification message and the data set metadata file.
- `"destinationS3Prefix"`: (Optional) The desired S3 prefix for the published data set,
  similar to a directory path in standard file systems. For example, if given the bucket name
  \"mybucket\" and the prefix \"myprefix/mydatasets\", the output file \"outputfile\" would
  be published to \"s3://mybucket/myprefix/mydatasets/outputfile\". If the prefix directory
  structure does not exist, it will be created. If no prefix is provided, the data set will
  be published to the S3 bucket root.
"""
function start_support_data_export(
    dataSetType,
    destinationS3BucketName,
    fromDate,
    roleNameArn,
    snsTopicArn;
    aws_config::AbstractAWSConfig=global_aws_config(),
)
    return marketplace_commerce_analytics(
        "StartSupportDataExport",
        Dict{String,Any}(
            "dataSetType" => dataSetType,
            "destinationS3BucketName" => destinationS3BucketName,
            "fromDate" => fromDate,
            "roleNameArn" => roleNameArn,
            "snsTopicArn" => snsTopicArn,
        );
        aws_config=aws_config,
        feature_set=SERVICE_FEATURE_SET,
    )
end
function start_support_data_export(
    dataSetType,
    destinationS3BucketName,
    fromDate,
    roleNameArn,
    snsTopicArn,
    params::AbstractDict{String};
    aws_config::AbstractAWSConfig=global_aws_config(),
)
    return marketplace_commerce_analytics(
        "StartSupportDataExport",
        Dict{String,Any}(
            mergewith(
                _merge,
                Dict{String,Any}(
                    "dataSetType" => dataSetType,
                    "destinationS3BucketName" => destinationS3BucketName,
                    "fromDate" => fromDate,
                    "roleNameArn" => roleNameArn,
                    "snsTopicArn" => snsTopicArn,
                ),
                params,
            ),
        );
        aws_config=aws_config,
        feature_set=SERVICE_FEATURE_SET,
    )
end


text:
"""
A `QuietRealTimePresenter` prints nothing as scenarios are being executed.
This is useful if you're only after a final report, after all features have been executed.
"""
struct QuietRealTimePresenter <: RealTimePresenter end
present(::QuietRealTimePresenter, ::Scenario) = nothing
present(::QuietRealTimePresenter, ::Gherkin.ScenarioStep) = nothing
present(::QuietRealTimePresenter, ::Gherkin.ScenarioStep, ::StepExecutionResult) = nothing
present(::QuietRealTimePresenter, ::Gherkin.Feature) = nothing

"""
This presenter prints a line to the console before, during, and after a scenario. Results are color
coded.
"""
struct ColorConsolePresenter <: RealTimePresenter
    io::IO
    colors::Dict{Type, Symbol}

    function ColorConsolePresenter(io::IO = stdout)
        colors = Dict{Type, Symbol}(
            NoStepDefinitionFound => :yellow,
            NonUniqueMatch => :magenta,
            SuccessfulStepExecution => :green,
            StepFailed => :red,
            UnexpectedStepError => :light_magenta,
            SkippedStep => :light_black
        )

        new(io, colors)
    end
end

"""
    stepformat(step::ScenarioStep)

Format a step according to Gherkin syntax.
"""
stepformat(step::Given) = "Given $(step.text)"
stepformat(step::When) = " When $(step.text)"
stepformat(step::Then) = " Then $(step.text)"

"""
    stepcolor(::Presenter, ::StepExecutionResult)

Get a console color for a given result when executing a scenario step.
"""
stepcolor(presenter::Presenter, step::StepExecutionResult) = presenter.colors[typeof(step)]

"""
    stepresultmessage(::StepExecutionResult)

A human readable message for a given result.
"""
stepresultmessage(step::StepFailed) = ["FAILED: " * step.assertion]
stepresultmessage(nomatch::NoStepDefinitionFound) = ["No match for '$(stepformat(nomatch.step))'"]
stepresultmessage(nonunique::NonUniqueMatch) = vcat(["Multiple matches found:"], ["  In " * location.filename for location in nonunique.locations])
stepresultmessage(::SuccessfulStepExecution) = []
stepresultmessage(::SkippedStep) = []
stepresultmessage(unexpected::UnexpectedStepError) = vcat(["Exception: $(string(unexpected.ex))"], ["  " * string(x) for x in unexpected.stack])
stepresultmessage(::StepExecutionResult) = []

function present(presenter::ColorConsolePresenter, feature::Feature)
    println()
    printstyled(presenter.io, "Feature: $(feature.header.description)\n"; color=:white)
end

function present(presenter::ColorConsolePresenter, scenario::Scenario)
    println()
    printstyled(presenter.io, "  Scenario: $(scenario.description)\n"; color=:blue)
end

function present(presenter::ColorConsolePresenter, step::Gherkin.ScenarioStep)
    printstyled(presenter.io, "    $(stepformat(step))"; color=:light_cyan)
end

function present(presenter::ColorConsolePresenter, step::Gherkin.ScenarioStep, result::StepExecutionResult)
    color = stepcolor(presenter, result)
    printstyled(presenter.io, "\r    $(stepformat(step))\n"; color=color)

    resultmessage = stepresultmessage(result)
    if !isempty(resultmessage)
        # The result message may be one or more lines. Indent them all.
        s = [" " ^ 8 * x * "\n" for x in resultmessage]
        println()
        println(join(s))
    end
end


text:
using TensorOperations
using KrylovKit
using LinearAlgebra

function one_site_dmrg( MPS::Dict,
                        MPO::Dict,
                        chi::Int64,
                        sweeps::Int64,
                        N::Int64)

    d = 2

    #MPS=init_random_MPS(d,chi,N)
    MPS=Normalize(MPS)

    L = construct_L(MPS, MPO, MPS, N)
    R = construct_R(MPS, MPO, MPS, N)

    global MPS,L,R

    for sweep in 1:Int(sweeps/2)
        for i in 1:N-1
            println(i)
            Energy,MPS[i],MPS[i+1] = optimize_one_site( MPS[i], MPS[i+1], MPO[i], L[i], R[i], "right", chi)
            println("Sweep ",sweep*2," Sites ",i," ",i+1," Energy ",Energy)
            L[i+1] = contract_from_left(L[i], MPS[i], MPO[i])
        end

        for i in N:-1:2
            Energy,MPS[i-1],MPS[i] = optimize_one_site( MPS[i-1], MPS[i], MPO[i], L[i], R[i], "left", chi)
            println("Sweep ",sweep*2,"Sites ",i," ",i+1,"Energy %-16.15f\n",Energy)
            R[i-1] = contract_from_right(R[i], MPS[i], MPO[i])
        end

    end

    return MPS

end

function optimize_one_site(A, B, W, E, F, dir, chi )
    function H_lin(v)
        ##     +--A--+
        ##     |  |  |
        ##     L--W--R
        ##     |  |  |
        #println("E ",size(E))
        #println("MPS ",size(v))
        #println("F ",size(F))
        #println("W ",size(W))
        @tensor temp_1[ :] := E[ 1, -2, -3] *  v[ -1, 1, -4]
        @tensor temp_2[:] :=  temp_1[-1,-2,-3,1] * F[1,-4,-5]
        @tensor temp[:] := temp_2[1, 2, -2, 3,-3] * W[ 2 3 ; 1 -1 ]
        return temp
    end

    if (dir == "right" )
        sA=size(A)

        en,V,info = eigsolve( H_lin , A ,  1  , :SR ; issymmetric = true  )
        en = en[1]
        V = V[1]
        sV = size(V)

        V=reshape(V,(sV[1]*sV[2],sV[3]))#s i j -> i*s, j
        A,S,V = svd(V,full=false)
        V=V'
        A = reshape( A , ( sA[1], sA[2], :) )
        #"ij,jl,slk->sik"
        S=diagm(S)
        @tensor B[:] := S[-2,1 ] * V[ 1,2 ] * B[-1,2,-3]
    else
        sB=size(B)
        print(sB)
        en,V,info = eigsolve( H_lin , B ,  1  , :SR ; issymmetric = true  )
        en = en[1]
        V = V[1]

        V=permutedims(V,(2,1,3))
        sV=size(V)
        V=reshape(V,(sV[1],sV[2]*sV[3]))  #s i j -> i s j -> i, s*j

        U,S,B = svd(V,full=false)
        B=B'
        B = permutedims(reshape(B,(:,sB[1],sB[3])),(2,1,3))
        #"sij,jk,kl->sil"
        S=diagm(S)
        @tensor A[:] := A[-1,-2,2] * U[ 2,3 ] * S[ 3,-3 ]
    end

    return en, A, B

end

max_stars_repo_name: aanniikkeettr/Tensor-Networks

max_issues_repo_name: aanniikkeettr/Tensor-Networks

max_forks_repo_name: aanniikkeettr/Tensor-Networks


text:
#=
Install script for the required packages
=#
using Pkg

println("Installing necessary packages. This may take some time.")
Pkg.add("DifferentialEquations")
Pkg.add("DynamicalSystems")
Pkg.add("ChaosTools")
Pkg.add("Plots")
Pkg.add("Pyplot")
Pkg.add("PlotlyJS")
Pkg.add("Images")
Pkg.add("ImageMagick")
Pkg.add("StatsBase")

max_stars_repo_name: csimal/chaos-and-fractals

max_issues_repo_name: csimal/chaos-and-fractals

max_forks_repo_name: csimal/chaos-and-fractals


text:
@testset "Strictly{LexicographicallyGreaterThan2Indicator}: $(fct_type), $(x_size) × $(y_size), $(T)" for fct_type in ["vector of variables", "vector affine function"], x_size in [2, 4], y_size in [2, 4], T in [Int, Float64]
    mock = MOIU.MockOptimizer(IndicatorMILPModel{T}())
    model = COIB.StrictlyLexicographicallyGreaterThan2Indicator{T}(mock)

    @test MOI.supports_constraint(model, MOI.VariableIndex, MOI.ZeroOne)
    @test MOI.supports_constraint(
        model,
        MOI.VectorAffineFunction{T},
        MOI.Indicator{MOI.ACTIVATE_ON_ONE, MOI.EqualTo{T}},
    )
    @test MOIB.supports_bridging_constraint(
        model,
        MOI.VectorAffineFunction{T},
        CP.Strictly{CP.LexicographicallyGreaterThan, T},
    )

    if T == Int
        x_array, _ = MOI.add_constrained_variables(model, [MOI.Integer() for _ in 1:(x_size * y_size)])
    elseif T == Float64
        x_array = MOI.add_variables(model, x_size * y_size)
    end
    x_matrix = reshape(x_array, x_size, y_size)

    fct = if fct_type == "vector of variables"
        MOI.VectorOfVariables(x_array)
    elseif fct_type == "vector affine function"
        MOIU.vectorize(x_array)
    else
        @assert false
    end
    c = MOI.add_constraint(model, fct, CP.Strictly{CP.LexicographicallyGreaterThan, T}(CP.LexicographicallyGreaterThan(y_size, x_size)))

    for i in 1:(x_size * y_size)
        @test MOI.is_valid(model, x_array[i])
    end
    @test MOI.is_valid(model, c)

    bridge = MOIBC.bridges(model)[MOI.ConstraintIndex{MOI.VectorOfVariables, CP.LexicographicallyGreaterThan}(-1)]

    @testset "Bridge properties" begin
        @test MOIBC.concrete_bridge_type(typeof(bridge), MOI.VectorOfVariables, CP.LexicographicallyGreaterThan) == typeof(bridge)
        @test MOIB.added_constrained_variable_types(typeof(bridge)) == [(MOI.ZeroOne,)]
        @test MOIB.added_constraint_types(typeof(bridge)) == [
            (MOI.ScalarAffineFunction{T}, MOI.EqualTo{T}),
            (MOI.VectorAffineFunction{T}, MOI.Indicator{MOI.ACTIVATE_ON_ONE, MOI.EqualTo{T}}),
            (MOI.VectorAffineFunction{T}, MOI.Indicator{MOI.ACTIVATE_ON_ONE, MOI.GreaterThan{T}}),
        ]

        @test MOI.get(bridge, MOI.NumberOfVariables()) == 2 * (x_size - 1) * y_size
        @test MOI.get(bridge, MOI.NumberOfConstraints{MOI.VariableIndex, MOI.ZeroOne}()) == 2 * (x_size - 1) * y_size
        @test MOI.get(bridge, MOI.NumberOfConstraints{MOI.ScalarAffineFunction{T}, MOI.EqualTo{T}}()) == x_size - 1 + (x_size - 1) * (y_size - 1)
        @test MOI.get(bridge, MOI.NumberOfConstraints{MOI.VectorAffineFunction{T}, MOI.Indicator{MOI.ACTIVATE_ON_ONE, MOI.EqualTo{T}}}()) == (x_size - 1) * y_size
        @test MOI.get(bridge, MOI.NumberOfConstraints{MOI.VectorAffineFunction{T}, MOI.Indicator{MOI.ACTIVATE_ON_ONE, MOI.GreaterThan{T}}}()) == (x_size - 1) * y_size

        @test Set(MOI.get(bridge, MOI.ListOfVariableIndices())) == Set([bridge.vars_eq..., bridge.vars_gt...])
        @test Set(MOI.get(bridge, MOI.ListOfConstraintIndices{MOI.VariableIndex, MOI.ZeroOne}())) == Set([bridge.vars_eq_bin..., bridge.vars_gt_bin...])
        @test Set(MOI.get(bridge, MOI.ListOfConstraintIndices{MOI.ScalarAffineFunction{T}, MOI.EqualTo{T}}())) == Set([vec(bridge.cons_move)..., bridge.cons_one_gt...])
        @test Set(MOI.get(bridge, MOI.ListOfConstraintIndices{MOI.VectorAffineFunction{T}, MOI.Indicator{MOI.ACTIVATE_ON_ONE, MOI.EqualTo{T}}}())) == Set(bridge.cons_indic_eq)
        @test Set(MOI.get(bridge, MOI.ListOfConstraintIndices{MOI.VectorAffineFunction{T}, MOI.Indicator{MOI.ACTIVATE_ON_ONE, MOI.GreaterThan{T}}}())) == Set(bridge.cons_indic_gt)
    end

    @testset "Sets of variables" begin
        @test length(bridge.vars_eq) == (x_size - 1) * y_size
        @test length(bridge.vars_gt) == (x_size - 1) * y_size

        for i in 1:(x_size - 1)
            # For each pair of columns (i, i+1)...
            @test length(bridge.vars_eq[i, :]) == y_size
            @test length(bridge.vars_gt[i, :]) == y_size

            for j in 1:y_size
                # For each pair of elements in the columns...
                @test MOI.is_valid(model, bridge.vars_eq[i, j])
                @test MOI.is_valid(model, bridge.vars_gt[i, j])
                @test MOI.is_valid(model, bridge.vars_eq_bin[i, j])
                @test MOI.is_valid(model, bridge.vars_gt_bin[i, j])

                @test MOI.get(model, MOI.ConstraintSet(), bridge.vars_eq_bin[i, j]) == MOI.ZeroOne()
                @test MOI.get(model, MOI.ConstraintSet(), bridge.vars_gt_bin[i, j]) == MOI.ZeroOne()
                @test MOI.get(model, MOI.ConstraintFunction(), bridge.vars_eq_bin[i, j]) == bridge.vars_eq[i, j]
                @test MOI.get(model, MOI.ConstraintFunction(), bridge.vars_gt_bin[i, j]) == bridge.vars_gt[i, j]
            end
        end
    end

    @testset "Constraint: exact one less-than pair of values per pair of columns" begin
        @test length(bridge.cons_one_gt) == (x_size - 1)

        for i in 1:(x_size - 1)
            # For each pair of columns (i, i+1)...
            @test MOI.is_valid(model, bridge.cons_one_gt[i])
            @test MOI.get(model, MOI.ConstraintSet(), bridge.cons_one_gt[i]) == MOI.EqualTo(one(T))

            f = MOI.get(model, MOI.ConstraintFunction(), bridge.cons_one_gt[i])
            @test length(f.terms) == y_size
            @test f.constant == zero(T)

            for j in 1:y_size
                t = f.terms[j]
                @test t.coefficient === one(T)
                @test t.variable == bridge.vars_gt[i, j]
            end
        end
    end

    @testset "Constraint: recursive definition starting at the end of the arrays" begin
        @test length(bridge.cons_move) == (x_size - 1) * (y_size - 1)

        for i in 1:(x_size - 1)
            # For each pair of columns (i, i+1)...
            @test length(bridge.cons_move[i, :]) == y_size - 1

            for j in 1:(y_size - 1)
                # For each pair of elements in the columns...
                @test MOI.is_valid(model, bridge.cons_move[i, j])
                @test MOI.get(model, MOI.ConstraintSet(), bridge.cons_move[i, j]) == MOI.EqualTo(zero(T))

                f = MOI.get(model, MOI.ConstraintFunction(), bridge.cons_move[i, j])
                @test length(f.terms) == 3
                @test f.constant == zero(T)

                t1 = f.terms[1]
                @test t1.coefficient === one(T)
                @test t1.variable == bridge.vars_eq[i, j]

                t2 = f.terms[2]
                @test t2.coefficient === -one(T)
                @test t2.variable == bridge.vars_eq[i, j + 1]

                t3 = f.terms[3]
                @test t3.coefficient === -one(T)
                @test t3.variable == bridge.vars_gt[i, j + 1]
            end
        end
    end

    @testset "Constraint: indicator for equality" begin
        @test length(bridge.cons_indic_eq) == (x_size - 1) * y_size

        for i in 1:(x_size - 1)
            # For each pair of columns (i, i+1)...
            @test length(bridge.cons_indic_eq[i, :]) == y_size

            for j in 1:y_size
                # For each pair of elements in the columns...
                @test MOI.is_valid(model, bridge.cons_indic_eq[i, j])
                @test MOI.get(model, MOI.ConstraintSet(), bridge.cons_indic_eq[i, j]) == MOI.Indicator{MOI.ACTIVATE_ON_ONE}(MOI.EqualTo(zero(T)))

                f = MOI.get(model, MOI.ConstraintFunction(), bridge.cons_indic_eq[i, j])
                @test length(f.terms) == 3
                @test f.constants == zeros(T, 2)

                t1 = f.terms[1]
                @test t1.output_index == 1
                @test t1.scalar_term.coefficient === one(T)
                @test t1.scalar_term.variable == bridge.vars_eq[i, j]

                t2 = f.terms[2]
                @test t2.output_index == 2
                @test t2.scalar_term.coefficient === one(T)
                @test t2.scalar_term.variable == x_matrix[i, j]

                t3 = f.terms[3]
                @test t3.output_index == 2
                @test t3.scalar_term.coefficient === -one(T)
                @test t3.scalar_term.variable == x_matrix[i + 1, j]
            end
        end
    end

    @testset "Constraint: indicator for inequality" begin
        @test length(bridge.cons_indic_gt) == (x_size - 1) * y_size

        for i in 1:(x_size - 1)
            # For each pair of columns (i, i+1)...
            @test length(bridge.cons_indic_gt[i, :]) == y_size

            for j in 1:y_size
                # For each pair of elements in the columns...
                @test MOI.is_valid(model, bridge.cons_indic_gt[i, j])
                @test MOI.get(model, MOI.ConstraintSet(), bridge.cons_indic_gt[i, j]) == MOI.Indicator{MOI.ACTIVATE_ON_ONE}(MOI.GreaterThan(zero(T)))

                f = MOI.get(model, MOI.ConstraintFunction(), bridge.cons_indic_gt[i, j])
                @test length(f.terms) == 3
                @test f.constants == zeros(T, 2)

                t1 = f.terms[1]
                @test t1.output_index == 1
                @test t1.scalar_term.coefficient === one(T)
                @test t1.scalar_term.variable == bridge.vars_gt[i, j]

                t2 = f.terms[2]
                @test t2.output_index == 2
                @test t2.scalar_term.coefficient === one(T)
                @test t2.scalar_term.variable == x_matrix[i, j]

                t3 = f.terms[3]
                @test t3.output_index == 2
                @test t3.scalar_term.coefficient === -one(T)
                @test t3.scalar_term.variable == x_matrix[i + 1, j]
            end
        end
    end
end


text:
module DaprProto

include("generated/appcallback/AppCallbackClients.jl")
include("generated/dapr/DaprClients.jl")
include("generated/operator/OperatorClients.jl")
include("generated/placement/PlacementClients.jl")
include("generated/sentry/SentryClients.jl")
include("generated/service_invocation/ServiceInvocationClients.jl")

end # module


text:
# Tests for deprecated functionality.
#
# These can't be run with --depwarn=error, so currently require special
# treatment when run inside the test system.

using Test
using Logging

module DeprecationTests # to test @deprecate
    f() = true

    # test the Symbol path of @deprecate
    @deprecate f1 f
    @deprecate f2 f false # test that f2 is not exported

    # test the Expr path of @deprecate
    @deprecate f3() f()
    @deprecate f4() f() false # test that f4 is not exported
    @deprecate f5(x::T) where T f()

    # test deprecation of a constructor
    struct A{T} end
    @deprecate A{T}(x::S) where {T, S} f()

    # test that @deprecate_moved can be overridden by an import
    Base.@deprecate_moved foo1234 "Foo"
    Base.@deprecate_moved bar "Bar" false
end # module
module Foo1234
    export foo1234
    foo1234(x) = x+1
end

# issue #21972
struct T21972
    @noinline function T21972()
        Base.depwarn("something", :T21972)
        new()
    end
end

@testset "@deprecate" begin
    using .DeprecationTests
    using .Foo1234
    @test foo1234(3) == 4
    @test_throws ErrorException DeprecationTests.bar(3)

    # 22845
    ex = :(module M22845; import ..DeprecationTests: bar;
                          bar(x::Number) = x + 3; end)
    @test_warn "importing deprecated binding" eval(ex)
    @test @test_nowarn(DeprecationTests.bar(4)) == 7

    # enable when issue #22043 is fixed
    # @test @test_warn "f1 is deprecated, use f instead." f1()
    # @test @test_nowarn f1()

    # @test_throws UndefVarError f2() # not exported
    # @test @test_warn "f2 is deprecated, use f instead." DeprecationTests.f2()
    # @test @test_nowarn DeprecationTests.f2()

    # @test @test_warn "f3() is deprecated, use f() instead." f3()
    # @test @test_nowarn f3()

    # @test_throws UndefVarError f4() # not exported
    # @test @test_warn "f4() is deprecated, use f() instead." DeprecationTests.f4()
    # @test @test_nowarn DeprecationTests.f4()

    # @test @test_warn "f5(x::T) where T is deprecated, use f() instead." f5(1)
    # @test @test_nowarn f5(1)

    # @test @test_warn "A{T}(x::S) where {T, S} is deprecated, use f() instead." A{Int}(1.)
    # @test @test_nowarn A{Int}(1.)

    # issue #21972
    @noinline function f21972()
        T21972()
    end
    @test_deprecated "something" f21972()
end

f24658() = depwarn24658()

depwarn24658() = Base.firstcaller(backtrace(), :_func_not_found_)

@testset "firstcaller" begin
    # issue #24658
    @test eval(:(if true; f24658(); end)) == (Ptr{Cvoid}(0),StackTraces.UNKNOWN)
end

# issue #25130
f25130() = Base.depwarn("f25130 message", :f25130)
# The following test is for the depwarn behavior of expressions evaluated at
# top-level, so we can't use the usual `collect_test_logs()` / `with_logger()`
testlogger = Test.TestLogger()
prev_logger = global_logger(testlogger)
# Each call at top level should be distinct. This won't be true if they're
# attributed to internal C frames (including generic dispatch machinery)
f25130()
f25130()
testlogs = testlogger.logs
@test length(testlogs) == 2
@test testlogs[1].id != testlogs[2].id
@test testlogs[1].kwargs.caller.func == Symbol("top-level scope")
@test all(l.message == "f25130 message" for l in testlogs)
global_logger(prev_logger)

#-------------------------------------------------------------------------------
# BEGIN 0.7 deprecations

module LogTest
    function bar(io)
        info(io,"barinfo")
        warn(io,"barwarn")
        Base.display_error(io,"barerror",backtrace())
    end
    function pooh(io)
        info(io,"poohinfo")
        warn(io,"poohwarn")
        Base.display_error(io,"pooherror",backtrace())
    end
end
function foo(io)
    info(io,"fooinfo")
    warn(io,"foowarn")
    Base.display_error(io,"fooerror",backtrace())
end

# Silence the flurry of depwarns for now.
with_logger(NullLogger()) do

@testset "Deprecated logging" begin

# Test info
@test contains(sprint(info, "test"), "INFO:")
@test contains(sprint(info, "test"), "INFO: test")
@test contains(sprint(info, "test ", 1, 2, 3), "INFO: test 123")
@test contains(sprint(io->info(io,"test", prefix="MYINFO: ")), "MYINFO: test")

# Test warn
@test contains(sprint(Base.warn_once, "test"), "WARNING: test")
@test isempty(sprint(Base.warn_once, "test"))

@test contains(sprint(warn), "WARNING:")
@test contains(sprint(warn, "test"), "WARNING: test")
@test contains(sprint(warn, "test ", 1, 2, 3), "WARNING: test 123")
@test contains(sprint(io->warn(io, "test", prefix="MYWARNING: ")), "MYWARNING: test")
@test contains(sprint(io->warn(io, "testonce", once=true)), "WARNING: testonce")
@test isempty(sprint(io->warn(io, "testonce", once=true)))
@test !isempty(sprint(io->warn(io, "testonce", once=true, key=hash("testonce",hash("testanother")))))
let bt = backtrace()
    ws = split(chomp(sprint(io->warn(io, "test", bt = bt))), '\n')
    bs = split(chomp(sprint(Base.show_backtrace, bt)), '\n')
    @test contains(ws[1],"WARNING: test")
    for (l,b) in zip(ws[2:end],bs[2:end])
        @test contains(l, b)
    end
end

# PR #16213
@test all(contains.(sprint(LogTest.bar), ["INFO: barinfo", "WARNING: barwarn", "ERROR: \"barerror\""]))
@test all(contains.(sprint(LogTest.pooh), ["INFO: poohinfo", "WARNING: poohwarn", "ERROR: \"pooherror\""]))
@test all(contains.(sprint(foo), ["INFO: fooinfo", "WARNING: foowarn", "ERROR: \"fooerror\""]))

logging(DevNull, LogTest, :bar;  kind=:info)
@test all(contains.(sprint(LogTest.bar), ["WARNING: barwarn", "ERROR: \"barerror\""]))
@test all(contains.(sprint(LogTest.pooh), ["INFO: poohinfo", "WARNING: poohwarn", "ERROR: \"pooherror\""]))
@test all(contains.(sprint(foo), ["INFO: fooinfo", "WARNING: foowarn", "ERROR: \"fooerror\""]))

logging(DevNull, LogTest;  kind=:info)
@test all(contains.(sprint(LogTest.bar), ["WARNING: barwarn", "ERROR: \"barerror\""]))
@test all(contains.(sprint(LogTest.pooh), ["WARNING: poohwarn", "ERROR: \"pooherror\""]))
@test all(contains.(sprint(foo), ["INFO: fooinfo", "WARNING: foowarn", "ERROR: \"fooerror\""]))

logging(DevNull;  kind=:info)
@test all(contains.(sprint(LogTest.bar), ["WARNING: barwarn", "ERROR: \"barerror\""]))
@test all(contains.(sprint(LogTest.pooh), ["WARNING: poohwarn", "ERROR: \"pooherror\""]))
@test all(contains.(sprint(foo), ["WARNING: foowarn", "ERROR: \"fooerror\""]))

logging(kind=:info)
@test all(contains.(sprint(LogTest.bar), ["INFO: barinfo", "WARNING: barwarn", "ERROR: \"barerror\""]))
@test all(contains.(sprint(LogTest.pooh), ["INFO: poohinfo", "WARNING: poohwarn", "ERROR: \"pooherror\""]))
@test all(contains.(sprint(foo), ["INFO: fooinfo", "WARNING: foowarn", "ERROR: \"fooerror\""]))

logging(DevNull, LogTest, :bar;  kind=:warn)
@test all(contains.(sprint(LogTest.bar), ["INFO: barinfo", "ERROR: \"barerror\""]))
@test all(contains.(sprint(LogTest.pooh), ["INFO: poohinfo", "WARNING: poohwarn", "ERROR: \"pooherror\""]))
@test all(contains.(sprint(foo), ["INFO: fooinfo", "WARNING: foowarn", "ERROR: \"fooerror\""]))

logging(DevNull, LogTest;  kind=:warn)
@test all(contains.(sprint(LogTest.bar), ["INFO: barinfo", "ERROR: \"barerror\""]))
@test all(contains.(sprint(LogTest.pooh), ["INFO: poohinfo", "ERROR: \"pooherror\""]))
@test all(contains.(sprint(foo), ["INFO: fooinfo", "WARNING: foowarn", "ERROR: \"fooerror\""]))

logging(DevNull;  kind=:warn)
@test all(contains.(sprint(LogTest.bar), ["INFO: barinfo", "ERROR: \"barerror\""]))
@test all(contains.(sprint(LogTest.pooh), ["INFO: poohinfo", "ERROR: \"pooherror\""]))
@test all(contains.(sprint(foo), ["INFO: fooinfo", "ERROR: \"fooerror\""]))

logging(kind=:warn)
@test all(contains.(sprint(LogTest.bar), ["INFO: barinfo", "WARNING: barwarn", "ERROR: \"barerror\""]))
@test all(contains.(sprint(LogTest.pooh), ["INFO: poohinfo", "WARNING: poohwarn", "ERROR: \"pooherror\""]))
@test all(contains.(sprint(foo), ["INFO: fooinfo", "WARNING: foowarn", "ERROR: \"fooerror\""]))

logging(DevNull, LogTest, :bar;  kind=:error)
@test all(contains.(sprint(LogTest.bar), ["INFO: barinfo", "WARNING: barwarn"]))
@test all(contains.(sprint(LogTest.pooh), ["INFO: poohinfo", "WARNING: poohwarn", "ERROR: \"pooherror\""]))
@test all(contains.(sprint(foo), ["INFO: fooinfo", "WARNING: foowarn", "ERROR: \"fooerror\""]))

logging(DevNull, LogTest;  kind=:error)
@test all(contains.(sprint(LogTest.bar), ["INFO: barinfo", "WARNING: barwarn"]))
@test all(contains.(sprint(LogTest.pooh), ["INFO: poohinfo", "WARNING: poohwarn"]))
@test all(contains.(sprint(foo), ["INFO: fooinfo", "WARNING: foowarn", "ERROR: \"fooerror\""]))

logging(DevNull;  kind=:error)
@test all(contains.(sprint(LogTest.bar), ["INFO: barinfo", "WARNING: barwarn"]))
@test all(contains.(sprint(LogTest.pooh), ["INFO: poohinfo", "WARNING: poohwarn"]))
@test all(contains.(sprint(foo), ["INFO: fooinfo", "WARNING: foowarn"]))

logging(kind=:error)
@test all(contains.(sprint(LogTest.bar), ["INFO: barinfo", "WARNING: barwarn", "ERROR: \"barerror\""]))
@test all(contains.(sprint(LogTest.pooh), ["INFO: poohinfo", "WARNING: poohwarn", "ERROR: \"pooherror\""]))
@test all(contains.(sprint(foo), ["INFO: fooinfo", "WARNING: foowarn", "ERROR: \"fooerror\""]))

logging(DevNull, LogTest, :bar)
@test sprint(LogTest.bar) == ""
@test all(contains.(sprint(LogTest.pooh), ["INFO: poohinfo", "WARNING: poohwarn", "ERROR: \"pooherror\""]))
@test all(contains.(sprint(foo), ["INFO: fooinfo", "WARNING: foowarn", "ERROR: \"fooerror\""]))

logging(DevNull, LogTest)
@test sprint(LogTest.bar) == ""
@test sprint(LogTest.pooh) == ""
@test all(contains.(sprint(foo), ["INFO: fooinfo", "WARNING: foowarn", "ERROR: \"fooerror\""]))

logging(DevNull)
@test sprint(LogTest.bar) == ""
@test sprint(LogTest.pooh) == ""
@test sprint(foo) == ""

logging()
@test all(contains.(sprint(LogTest.bar), ["INFO: barinfo", "WARNING: barwarn", "ERROR: \"barerror\""]))
@test all(contains.(sprint(LogTest.pooh), ["INFO: poohinfo", "WARNING: poohwarn", "ERROR: \"pooherror\""]))
@test all(contains.(sprint(foo), ["INFO: fooinfo", "WARNING: foowarn", "ERROR: \"fooerror\""]))

end # @testset

end

# END 0.7 deprecations

max_stars_repo_name: YingboMa/julia

max_issues_repo_name: YingboMa/julia

max_forks_repo_name: YingboMa/julia


text:
export bubble_sort!, bubble_sort_perm

function bubble_sort_pass!(v::AbstractVector)
    issorted = true
    for i in 1:(length(v)-1)
        if v[i] > v[i+1]
            v[i],v[i+1]=v[i+1],v[i]
            issorted = false
        end
    end
    return issorted
end

function bubble_sort!(v::AbstractVector)
    issorted = false
    while !issorted
        issorted = bubble_sort_pass!(v)
    end
    v
end

function bubble_sort_pass!(v::AbstractVector, fixed::AbstractVector)
    length(fixed) == length(v) || error("v and fixed must be the same length")
    issorted = true
    for i in 1:(length(v)-1)
        x,nxt = i,i+1
        ready = false
        while !ready
            nxt > length(v) && break
            if v[x] > v[nxt]
                if fixed[x] | fixed[nxt]
                    nxt = nxt+=1
                    continue
                end
                v[x],v[nxt]  = v[nxt],v[x]
                issorted = false
                ready = true
            else
                ready = true
            end
        end
    end
    return issorted
end

function bubble_sort_perm_pass!(v::AbstractVector, fixed::AbstractVector, by::AbstractVector)
    length(fixed) == length(v) || error("v and fixed must be the same length")
    issorted = true
    for i in 1:(length(v)-1)
        x,nxt = i,i+1
        ready = false
        while !ready
            nxt > length(v) && break
            if by[x] > by[nxt]
                if fixed[x] | fixed[nxt]
                    nxt = fixed[nxt] ? nxt+=1 : nxt
                    x   = fixed[x]   ? x+=1   : x
                    continue
                end
                v[x],v[nxt]  = v[nxt],v[x]
                by[x],by[nxt]= by[nxt],by[x]
                issorted = false
                ready = true
            else
                ready = true
            end
        end
    end
    return issorted
end

function bubble_sort!(v::AbstractVector, fixed::AbstractVector)
    issorted = false
    while !issorted
        issorted = bubble_sort_pass!(v,fixed)
    end
    v
end

function bubble_sort_perm(v::AbstractVector, fixed::AbstractVector,perm = collect(1:length(v)))
    issorted = false
    v2 = copy(v)
    while !issorted
        issorted = bubble_sort_perm_pass!(perm,fixed,v2)
    end
    perm
end


text:
using ImageFiltering, ImageCore
using Test

@testset "gradient" begin
    y, x = 1:5, (1:7)'
    img_grads = ((y .* fill(1, size(x)), 1, 0),
            (fill(1, size(y)) .* x, 0, 1),
            (map(Gray, y.*fill(1.0, size(x))), Gray(1), Gray(0)),
            (map(v->RGB(v,0,0), y.*fill(1.0, size(x))), RGB(1,0,0), RGB(0,0,0)))
    for (img, ey, ex) in img_grads
        for kernelfunc in (KernelFactors.ando3, KernelFactors.sobel, KernelFactors.prewitt,
                           KernelFactors.ando4, KernelFactors.ando5, KernelFactors.bickley,
                           KernelFactors.scharr,
                           Kernel.ando3, Kernel.sobel, Kernel.prewitt,
                           Kernel.ando4, Kernel.ando5, Kernel.scharr, Kernel.bickley)
            gy, gx = @inferred(imgradients(img, kernelfunc, Inner()))
            for val in gy
                @test _abs(val - ey) < 1e-4
            end
            for val in gx
                @test _abs(val - ex) < 1e-4
            end
            gy, gx = imgradients(img, kernelfunc, Pad(:replicate))
            @test axes(gy) == axes(gx) == axes(img)
        end
        for funcpairs in ((Kernel.ando3, KernelFactors.ando3),
                          (Kernel.sobel, KernelFactors.sobel),
                          (Kernel.prewitt, KernelFactors.prewitt),
                          (Kernel.scharr, KernelFactors.scharr),
                          (Kernel.bickley, KernelFactors.bickley))
            fk, fkf = funcpairs
            ky, kx = fk()
            gmy, gmx = imfilter(img, ky), imfilter(img, kx)
            gy, gx = imgradients(img, fkf)
            @test isapprox(gmy, gy; atol=1e-8)
            @test isapprox(gmx, gx; atol=1e-8)
        end
        for fk in (Kernel.ando4, Kernel.ando5)
            ky, kx = fk()
            gy, gx = imfilter(img, ky, Inner()), imfilter(img, kx, Inner())
            for val in gy
                @test _abs(val - ey) < 1e-4
            end
            for val in gx
                @test _abs(val - ex) < 1e-4
            end
        end
    end
    # 3d
    y, x, z = 1:5, (1:7)', reshape(1:6, 1, 1, 6)
    for (img, ey, ex, ez) in ((y .* fill(1, size(x)) .* fill(1, size(z)), 1, 0, 0),
                              (fill(1, size(y)) .* x .* fill(1, size(z)), 0, 1, 0),
                              (fill(1, size(y)) .* fill(1, size(x)) .* z, 0, 0, 1))
        for kernelfunc in (KernelFactors.ando3, KernelFactors.sobel, KernelFactors.prewitt,
                           KernelFactors.scharr, KernelFactors.bickley,
                           Kernel.ando3, Kernel.sobel, Kernel.prewitt,
                           Kernel.scharr, Kernel.bickley)
            gy, gx, gz = @inferred(imgradients(img, kernelfunc, Inner()))
            for val in gy
                @test _abs(val - ey) < 1e-4
            end
            for val in gx
                @test _abs(val - ex) < 1e-4
            end
            for val in gz
                @test _abs(val - ez) < 1e-4
            end
        end
    end
end

nothing


text:
function score!(cur_player::AbstractPlayer, cur_label::AbstractString, cur_value::Int)

  # --------------------
  #  verify score comes
  #   from valid label
  # --------------------

  in(cur_label, allowed_labels) ||
    error("Invalid label in score call: $(cur_label).")

  # ------------------
  #  handle soft drop
  # ------------------

  if cur_label == "soft"
    cur_player.score += 1
    cur_player.stash += 1

    return 1
  end

  # ---------------------
  #  calculate the score
  # ---------------------

  cur_dict = OrderedDict(
    "level" => cur_player.level,
    "lines" => cur_player.lines,
    "combo" => cur_player.combo,
    "stash" => cur_player.stash,
    "glues" => cur_player.glues,
    "total" => cur_player.score,
    "chain" => cur_player.chain,
    "piece" => cur_player.piece.name
  )

  cur_function = getfield(
    Tetris,
    Symbol("_score_$(cur_label)")
  )

  cur_score = cur_function(cur_value, cur_dict)

  # -------------------
  #  log current score
  # -------------------

  cur_round = cur_player.round

  if cur_round.is_keeping_score
    cur_position, cur_round.log_count =
      [0, 1] .+ cur_round.log_count

    cur_log = Log(
      cur_position,
      cur_label,
      cur_value,
      cur_score,
      values(cur_dict)...
    )

    push!(cur_round.logs, cur_log)
  end

  # ---------------------
  #  update player score
  # ---------------------

  if cur_label == "glue" || cur_label == "hold"
    cur_player.stash = 0
  else
    cur_player.score += cur_score
  end

  if is_repl
    cur_string = []
    push!(cur_string, crayon_dict["grey"])
    push!(cur_string, "\x1b[$(14);$(11)H")
    push!(cur_string, lpad(cur_player.score, 8, "0"))
    push!(cur_string, inv(crayon_dict["grey"]))
    push!(cur_string, "\x1b[u")
    print(cur_string...)
  end

  return cur_score

end

function _score_clear(cleared_count::Int, cur_dict::OrderedDict)
  if cleared_count == 4 && cur_dict["chain"] == "tetris"
    cur_score = 1200
  else
    cur_score = score_array[cleared_count]
  end

  cur_score *= cur_dict["level"]

  cur_score
end

function _score_combo(cleared_count::Int, cur_dict::OrderedDict)
  cur_score = ( cleared_count == 1 ) ? 20 : 50

  cur_score *= cur_dict["combo"]

  cur_score *= cur_dict["level"]

  cur_score
end

function _score_hard(drop_count::Int, cur_dict::OrderedDict)
  cur_score = drop_count

  cur_score *= 2

  cur_score
end

function _score_glue(stub_value::Int, cur_dict::OrderedDict)
  cur_score = cur_dict["stash"]

  cur_score
end

function _score_hold(stub_value::Int, cur_dict::OrderedDict)
  cur_score = cur_dict["stash"]

  cur_score
end


text:
"""
    function destriper(N::Int, dataset_baselines; n_iter=10)
# Brief description
You have to use this function to clean a noisy-TOD (comphrensive
of white and one-over-f nosie).

Suggestion for raster scan:
In order to suppress the atmospheric noise you should use many
baselines, better if you have observed those maps with different
scanning strategies. Now you can use only Left-Right and Up-Down
scanning.

The descriper solve the noisy-tod linear system:
[Pᵀ N⁻² P] m = [Pᵀ N⁻²] d

We can reduce it at the tipical solution of
A x = b where:

A = [Pᵀ N⁻² P]
b = [Pᵀ N⁻²] d

m = [Pᵀ N⁻² P]⁻¹ [Pᵀ N⁻²] d
These matris are too big, so we can you the conjugate gradient method
in order to solve the system.
"""

function destriper(N::Int, dataset_baselines; tolerance=1e-6, n_iter=100)
    b = zeros(N*N)
    for i in dataset_baselines
            tod = denoise(i.time_order_data, i.noise)
            a, h = binned_map(N, tod, i.pointing)
            b .+= reshape(a, N*N, 1)[:, 1]
    end

    function A(x)
        res = zeros(N, N)
        for i in dataset_baselines
            tod = get_tod(reshape(x, N, N), i.pointing)
            tod = denoise(tod, i.noise)
            mat, hmat = binned_map(N, tod, i.pointing)
            res .+= mat
        end
        return reshape(res, N*N, 1)[:, 1]
    end

    return ConjGrad.cg(A, b, tol=1e-10, maxIter=100)
    #return conjgrad(A, b, tol=1e-10, maxiter=100)

end

"""
    function denoise(tod::Array{Float64, 1}, noise_s::Array{Float64, 1})

# Brief Description
This function apply the N⁻² variance matrix to the TOD working in the
fourier space.
"""
function denoise(tod::Array{Float64, 1}, noise_s::Array{Float64, 1})
        FFTW.set_num_threads(12)
        ftod = fft(tod)
        ftod ./= noise_s
        tod = real.(ifft(ftod))
        return tod
end

"""
    function conjgrad(A, b; x0=0, tol=1e-6, maxiter=100)
# Brief description
"""
genblas_dot(x, y) = dot(x,y)
genblas_scal!(a, x) = x .*= a
genblas_axpy!(a, x, y) = y .+= a.*x
genblas_nrm2(x) = norm(x)
function conjgrad(A, b; tol=1e-6, maxiter=100)

    x = b .* 0
    r = b

    n = length(b)
    z = copy(r)

    p = z
    err = Inf
    i = 0

    residual_0 = genblas_nrm2(r)

    for i = 1:maxiter
        Ap = A(p)
        gamma = genblas_dot(r,z)
        alpha = gamma / genblas_dot(p, Ap)

        if alpha == Inf || alpha < 0
            return -13, iter
        end

        genblas_axpy!(alpha, p, x)
        genblas_axpy!(-alpha, Ap, r)

        residual = genblas_nrm2(r) / residual_0

        if residual <= tol
            println("Resi: ",residual)
            println("Iter: ",i)
            return x
        end

        z = copy(r)

        beta = genblas_dot(z, r) / gamma
        genblas_scal!(beta, p)
        genblas_axpy!(1.0, z, p)

        #println("Err: ", i, " ", residual)
    end

    return x

end

function make_dataset(tod::Array{Float64, 1}, pix_idx::Array{Int64, 1},
                      Num_bas::Int64
                      )
    len_tod = length(tod)
    len2    = Int(round(len_tod / Num_bas))

    mat_tod = reshape(tod, (len2, Num_bas))
    mat_point = reshape(pix_idx, (len2, Num_bas))

    dataset_baselines = Array{Tpn, 1}(undef, Num_bas)

    noise_spec = tod_noise(len2,  0.0083, 0.1, 3.0, 20.76)
    # Ora e` tutto uguale per ogni baselines, l'ideale sarebbe
    # inserire un fit, con il modello scelto, per ogni baselines

    for i in 1:Num_bas
        dat = Tpn(mat_tod[:, i], mat_point[:, i], noise_spec)
        dataset_baselines[i] = dat
    end

    return dataset_baselines

end


text:
module DataPacket

using LibSerialPort

export open, write_data, retrieve_data
export SP_MODE_READ

function write_data(data::String, path_to_datafile::String)
	open(path_to_datafile, "a") do f
		write(f, data)
	end
end

function retrieve_data(sp::SerialPort)
	datapacket = ""
	while true
		datapacket *= String(nonblocking_read(sp))
		if length(datapacket) >= 26
			println("Got a data packet:")
			print(datapacket)
			break
		end
		sleep(0.001)
	end
	return datapacket
end

end

max_stars_repo_name: alflhomme/lorasensors

max_issues_repo_name: alflhomme/lorasensors

max_forks_repo_name: alflhomme/lorasensors


text:
#=

  Added corner puzzle in constraint

  Problem from http://www.delphiforfun.org/Programs/AddedCorners.htm
  """
  This puzzle requires that you enter the digits 1 through 8 in the circles and
  squares (one digit in each figure) so that the number in each square is equal
  to the sum on the numbers in the circles which  adjoin it.
  ...

    C F C
    F   F
    C F C

  """

  Model created by Hakan Kjellerstrand, hakank@gmail.com
  See also my Julia page: http://www.hakank.org/julia/

=#

using ConstraintSolver, JuMP
using Cbc, GLPK, Ipopt
using Printf
const CS = ConstraintSolver
include("constraints_utils.jl")

function added_corner(print_solutions=true,all_solutions=true,timeout=6)

    cbc_optimizer = optimizer_with_attributes(Cbc.Optimizer, "logLevel" => 0)
    glpk_optimizer = optimizer_with_attributes(GLPK.Optimizer)
    ipopt_optimizer = optimizer_with_attributes(Ipopt.Optimizer)

    model = Model(optimizer_with_attributes(CS.Optimizer,   "all_solutions"=> all_solutions,
                                                            # "all_optimal_solutions"=>all_solutions,
                                                            "logging"=>[],

                                                            "traverse_strategy"=>:BFS,
                                                            # "traverse_strategy"=>:DFS,
                                                            # "traverse_strategy"=>:DBFS,

                                                            # "branch_split"=>:Smallest,
                                                            # "branch_split"=>:Biggest,
                                                            "branch_split"=>:InHalf,

                                                            # https://wikunia.github.io/ConstraintSolver.jl/stable/options/#branch_strategy-(:Auto)
                                                            "branch_strategy" => :IMPS, # default
                                                            # "branch_strategy" => :ABS, # Activity Based Search
                                                            # "activity.decay" => 0.999, # default 0.999
                                                            # "activity.max_probes" => 10, # default, 10
                                                            # "activity.max_confidence_deviation" => 20, # default 20

                                                            # "simplify"=>false,
                                                            # "simplify"=>true, # default

                                                            "time_limit"=>timeout,

                                                            # "backtrack" => false, # default true
                                                            # "backtrack_sorting" => false, # default true

                                                            # "lp_optimizer" => cbc_optimizer,
                                                            # "lp_optimizer" => glpk_optimizer,
                                                            # "lp_optimizer" => ipopt_optimizer,
                                        ))

    n = 8
    @variable(model, 1 <= x[1:n] <= n, Int)
    A,B,C,D,E,F,G,H = x
    @constraint(model, x in CS.AllDifferent())

    @constraint(model, B == A + C)
    @constraint(model, D == A + F)
    @constraint(model, E == C + H)
    @constraint(model, G == F + H)

    # Solve the problem
    optimize!(model)

    status = JuMP.termination_status(model)
    # println("status:$status")
    num_sols = 0
    if status == MOI.OPTIMAL
        num_sols = MOI.get(model, MOI.ResultCount())
        println("num_sols:$num_sols\n")
        if print_solutions
            for sol in 1:num_sols
                println("solution #$sol")
                x_val = convert.(Integer,JuMP.value.(x; result=sol))
                println("x:$x_val")
                @printf "%d %d %d\n" x_val[1] x_val[2] x_val[3]
                @printf "%d   %d\n" x_val[4] x_val[5]
                @printf "%d %d %d\n" x_val[6] x_val[7] x_val[8]
                println()

            end
        end
    else
        println("status:$status")
    end

    return status, num_sols
end

@time added_corner()

max_stars_repo_name: tias/hakank

max_issues_repo_name: tias/hakank

max_forks_repo_name: tias/hakank


text:
module Experimental
using Setfield
using ConstructionBase: constructorof
import Setfield: get, set
export MultiPropertyLens

const NNamedTupleLens{N,s} = NamedTuple{s, T} where {T <: NTuple{N, Lens}}
struct MultiPropertyLens{L <: NNamedTupleLens} <: Lens
    lenses::L
end

_keys(::Type{MultiPropertyLens{NamedTuple{s,T}}}) where {s,T} = s
@generated function get(obj, l::MultiPropertyLens)
    get_arg(fieldname) = :($fieldname = get(obj.$fieldname, l.lenses.$fieldname))
    args = map(get_arg, _keys(l))
    Expr(:tuple, args...)
end

@generated function set(obj, l::MultiPropertyLens, val)
    T = obj
    args = map(fieldnames(T)) do fn
        if fn in _keys(l)
            quote
                obj_inner = obj.$fn
                lens_inner = l.lenses.$fn
                val_inner = val.$fn
                set(obj_inner, lens_inner, val_inner)
            end
        else
            :(obj.$fn)
        end
    end
    Expr(:block,
        Expr(:meta, :inline),
        Expr(:call, :(constructorof($T)), args...)
    )
end

function Base.show(io::IO, l::MultiPropertyLens)
    print(io, "MultiPropertyLens(")
    print(io, l.lenses)
    print(io, ')')
end
end


text:
const default_initializer = ArrayInitializer(t -> t<:Union{Tuple, NamedTuple, Pair}, (T, sz) -> similar(arrayof(T), sz))

"""
    collect_columns(itr)

Collect an iterable as a `Columns` object if it iterates `Tuples` or `NamedTuples`, as a normal
`Array` otherwise.

# Examples

    s = [(1,2), (3,4)]
    collect_columns(s)

    s2 = Iterators.filter(isodd, 1:8)
    collect_columns(s2)
"""
collect_columns(itr) = collect_structarray(itr, initializer = default_initializer)
collect_empty_columns(itr) = collect_empty_structarray(itr, initializer = default_initializer)

grow_to_columns!(args...) = grow_to_structarray!(args...)
collect_to_columns!(args...) = collect_to_structarray!(args...)

function collect_columns_flattened(itr)
    elem = iterate(itr)
    @assert elem !== nothing
    el, st = elem
    collect_columns_flattened(itr, el, st)
end

function collect_columns_flattened(itr, el, st)
    while isempty(el)
        elem = iterate(itr, st)
        elem === nothing && return collect_empty_columns(el)
        el, st = elem
    end
    dest = collect_columns(el)
    collect_columns_flattened!(dest, itr, el, st)
end

function collect_columns_flattened!(dest, itr, el, st)
    while true
        elem = iterate(itr, st)
        elem === nothing && break
        el, st = elem
        dest = grow_to_columns!(dest, el)
    end
    return dest
end

function collect_columns_flattened(itr, el::Pair, st)
    while isempty(el.second)
        elem = iterate(itr, st)
        elem === nothing && return collect_empty_columns(el.first => i for i in el.second)
        el, st = elem
    end
    dest_data = collect_columns(el.second)
    dest_key = collect_columns(el.first for i in dest_data)
    collect_columns_flattened!(Columns(dest_key => dest_data), itr, el, st)
end

function collect_columns_flattened!(dest::Columns{<:Pair}, itr, el::Pair, st)
    dest_key, dest_data = columns(dest)
    while true
        elem = iterate(itr, st)
        elem === nothing && break
        el, st = elem
        n = length(dest_data)
        dest_data = grow_to_columns!(dest_data, el.second)
        dest_key = grow_to_columns!(dest_key, el.first for i in (n+1):length(dest_data))
    end
    return Columns(dest_key => dest_data)
end


text:
module DS

export data_splitting, DS_CI

using Random
using Distributions
using IntervalSets
using LinearAlgebra: I
include("../data.jl")
using .DataUtil
include("sfs.jl")
using .SFS

function data_splitting(𝐗, 𝐲, k, 𝚺; ratio=0.5)
    # data splitting
    idx = shuffle(1:length(𝐲)) # randomize data
    𝐗_FS, 𝐗_SI = train_val_split(𝐗[idx, :]; ratio=ratio)
    𝐲_FS, 𝐲_SI = train_val_split(𝐲[idx]; ratio=ratio)

    # feature selection
    𝐴_obs, _ = sfs(𝐗_FS, 𝐲_FS, k)

    # statistical inference (hypothesis test)
    𝐗_obs = 𝐗_SI[:, 𝐴_obs]
    𝜼 = 𝐗_obs / (𝐗_obs' * 𝐗_obs)
    p_values = Vector{Float64}(undef, length(𝐴_obs)) # p-value list
    for j = 1:length(𝐴_obs)
        𝜼ⱼ = 𝜼[:, j]
        z_obs = 𝜼ⱼ' * 𝐲_SI
        σ²_z = 𝜼ⱼ' * 𝚺 * 𝜼ⱼ # variance of z_obs
        𝜋_j = cdf(Normal(0, √(σ²_z)), z_obs)
        p_values[j] = 2 * min(𝜋_j, 1 - 𝜋_j)
    end
    return (𝐴_obs, p_values)
end

function DS_CI(𝐗, 𝐲, k, 𝚺; ratio=0.5, α=0.05)
    # data splitting
    idx = shuffle(1:length(𝐲)) # randomize data
    𝐗_FS, 𝐗_SI = train_val_split(𝐗[idx, :]; ratio=ratio)
    𝐲_FS, 𝐲_SI = train_val_split(𝐲[idx]; ratio=ratio)

    # feature selection
    𝐴_obs, _ = sfs(𝐗_FS, 𝐲_FS, k)

    # statistical inference (estimate confidence interval)
    Φ = cquantile(Normal(0, 1), α / 2)
    𝐗_obs = 𝐗_SI[:, 𝐴_obs]
    𝜼 = 𝐗_obs / (𝐗_obs' * 𝐗_obs)
    CIs = Vector{ClosedInterval{Float64}}(undef, length(𝐴_obs)) # CI list
    for j = 1:length(𝐴_obs)
        𝜼ⱼ = 𝜼[:, j]
        z_obs = 𝜼ⱼ' * 𝐲_SI
        σ_z = √(𝜼ⱼ' * 𝚺 * 𝜼ⱼ) # standard deviation of z_obs
        CIs[j] = (z_obs - σ_z * Φ)..(z_obs + σ_z * Φ)
    end
    return (𝐴_obs, CIs)
end

function FPR_experiment(β₀, 𝛃, N, k; α=0.05, denominator_upper=100, iter_size=5, ratio=0.5)
    H_0 = filter(i -> iszero(𝛃[i]), 1:length(𝛃)) # {j ∈ {1,...,d} | H_{0,j} is true}
    FPR = Array{Float64}(undef, length(H_0), iter_size)
    for l = 1:iter_size
        numerator = zeros(Int, length(H_0))
        denominator = zeros(Int, length(H_0))
        for _ = 1:denominator_upper
            𝐗, 𝐲_obs = centering(make_dataset(β₀, 𝛃, N)...)
            𝐴_obs, p_values = data_splitting(𝐗, 𝐲_obs, k, I; ratio=ratio)
            for (j, p_value) in zip(𝐴_obs, p_values)
                for (i, j′) in enumerate(H_0)
                    if j == j′ # H_{0,j} is true
                        numerator[i] += p_value < α ? 1 : 0
                        denominator[i] += 1
                        break
                    end
                end
            end
        end
        FPR[:, l] = numerator ./ denominator
    end
    return FPR
end

function TPR_experiment(β₀, 𝛃, N, k; α=0.05, denominator_upper=100, iter_size=5, ratio=0.5)
    H_1 = filter(i -> !iszero(𝛃[i]), 1:length(𝛃)) # {j ∈ {1,...,d} | H_{1,j} is true}
    TPR = Array{Float64}(undef, length(H_1), iter_size)
    for l = 1:iter_size
        numerator = zeros(Int, length(H_1))
        denominator = zeros(Int, length(H_1))
        for _ = 1:denominator_upper
            𝐗, 𝐲_obs = centering(make_dataset(β₀, 𝛃, N)...)
            𝐴_obs, p_values = data_splitting(𝐗, 𝐲_obs, k, I; ratio=ratio)
            for (j, p_value) in zip(𝐴_obs, p_values)
                for (i, j′) in enumerate(H_1)
                    if j == j′ # H_{1,j} is true
                        numerator[i] += p_value < α ? 1 : 0
                        denominator[i] += 1
                        break
                    end
                end
            end
        end
        TPR[:, l] = numerator ./ denominator
    end
    return TPR
end

end # module DS


text:
using Dates
# using GitCommand
using GitHub
using JSON
using Pkg
using Printf
using RegistryCI
using Test
using TimeZones

const AutoMerge = RegistryCI.AutoMerge

include("automerge-integration-utils.jl")

AUTOMERGE_INTEGRATION_TEST_REPO = ENV["AUTOMERGE_INTEGRATION_TEST_REPO"]::String
TEST_USER_GITHUB_TOKEN = ENV["BCBI_TEST_USER_GITHUB_TOKEN"]::String
INTEGRATION_TEST_READ_ONLY_TOKEN = ENV["INTEGRATION_TEST_READ_ONLY_TOKEN"]::String
GIT = "git"
auth = GitHub.authenticate(TEST_USER_GITHUB_TOKEN)
whoami = RegistryCI.AutoMerge.username(GitHub.DEFAULT_API, auth)
repo_url_without_auth = "https://github.com/$(AUTOMERGE_INTEGRATION_TEST_REPO)"
repo_url_with_auth = "https://$(whoami):$(TEST_USER_GITHUB_TOKEN)@github.com/$(AUTOMERGE_INTEGRATION_TEST_REPO)"
repo = GitHub.repo(AUTOMERGE_INTEGRATION_TEST_REPO; auth = auth)
@test success(`$(GIT) --version`)
@info("Authenticated to GitHub as \"$(whoami)\"")

_delete_branches_older_than = Dates.Hour(3)
delete_old_pull_request_branches(
    repo_url_with_auth,
    _delete_branches_older_than;
    GIT = GIT,
)

requires_commit = "1c843ed4d4d568345ca557fea48f43efdcd0271f"
hello_world_commit1 = "197e0e03f3f840f830cb5095ab6407d00fbee61c"
hello_world_commit2 = "57b0aec49622faa962c6752d4bc39a62b91fe37c"

@testset "Integration tests" begin
    for (test_number, (master_dir, feature_dir, public_dir, title, check_license, pass, commit)) in enumerate([
            ("master_1", "feature_1", "", "New package: Requires v1.0.0", true, true, requires_commit), # OK: new package
            ("master_1", "feature_1", "", "New package: Requires v1.0.0", true, false, "659e09770ba9fda4a503f8bf281d446c9583ff3b"), # FAIL: wrong commit!
            ("master_2", "feature_2", "", "New version: Requires v2.0.0", false, true, requires_commit),            # OK: new version
            ("master_1", "feature_3", "", "New package: Req v1.0.0", false, false, requires_commit),                # FAIL: name too short
            ("master_2", "feature_4", "", "New version: Requires v2.0.1", false, false, requires_commit),           # FAIL: skips v2.0.0
            ("master_3", "feature_5", "", "New version: Requires v2.0.0", false, false, requires_commit),           # FAIL: modifies extra file
            ("master_1", "feature_6", "", "New package: HelloWorldC_jll v1.0.6+0", false, true, hello_world_commit1),   # OK: new JLL package
            ("master_4", "feature_7", "", "New version: HelloWorldC_jll v1.0.8+0", false, true, hello_world_commit2),   # OK: new JLL version
            ("master_1", "feature_8", "", "New package: HelloWorldC_jll v1.0.6+0", false, false, hello_world_commit1),  # FAIL: unallowed dependency
            ("master_1", "feature_1", "public_1", "New package: Requires v1.0.0", false, true, requires_commit),    # OK: no UUID conflict
            ("master_1", "feature_1", "public_2", "New package: Requires v1.0.0", false, false, requires_commit),  # FAIL: UUID conflict, name differs
            ("master_1", "feature_1", "public_3", "New package: Requires v1.0.0", false, false, requires_commit),  # FAIL: UUID conflict, repo differs
            ("master_1", "feature_1", "public_4", "New package: Requires v1.0.0", false, true, requires_commit),   # OK: UUID conflict but name and repo match
        ])
        @info "Performing integration tests with settings" test_number master_dir feature_dir title pass
        with_master_branch(templates(master_dir), "master"; GIT = GIT, repo_url = repo_url_with_auth) do master
            with_feature_branch(templates(feature_dir), master; GIT = GIT, repo_url = repo_url_with_auth) do feature
                public_registries = String[]
                if public_dir != ""
                    public_git_repo = generate_public_registry(public_dir, GIT)
                    push!(public_registries, "file://$(public_git_repo)/.git")
                end
                head = feature
                base = master
                body = """
                - Foo: Bar
                - Commit: $commit
                - Hello: World
                """
                params = Dict("title" => title,
                              "head" => head,
                              "base" => base,
                              "body" => body)

                sleep(1)
                pr = GitHub.create_pull_request(repo; auth = auth, params = params)
                pr = wait_pr_compute_mergeability(GitHub.DEFAULT_API, repo, pr; auth = auth)
                @test pr.mergeable
                sleep(1)
                with_pr_merge_commit(pr, repo_url_without_auth; GIT = GIT) do build_dir
                    withenv("AUTOMERGE_GITHUB_TOKEN" => TEST_USER_GITHUB_TOKEN,
                            "TRAVIS_BRANCH" => master,
                            "TRAVIS_BUILD_DIR" => build_dir,
                            "TRAVIS_EVENT_TYPE" => "pull_request",
                            "TRAVIS_PULL_REQUEST" => string(pr.number),
                            "TRAVIS_PULL_REQUEST_SHA" => string(AutoMerge.pull_request_head_sha(pr)),
                            "TRAVIS_REPO_SLUG" => AUTOMERGE_INTEGRATION_TEST_REPO) do
                        sleep(1)
                        run_thunk = () -> AutoMerge.run(;
                                      merge_new_packages = true,
                                      merge_new_versions = true,
                                      new_package_waiting_period = Minute(typemax(Int32)),
                                      new_jll_package_waiting_period = Minute(typemax(Int32)),
                                      new_version_waiting_period = Minute(typemax(Int32)),
                                      new_jll_version_waiting_period = Minute(typemax(Int32)),
                                      registry = AUTOMERGE_INTEGRATION_TEST_REPO,
                                      authorized_authors = String[whoami],
                                      authorized_authors_special_jll_exceptions = String[whoami],
                                      error_exit_if_automerge_not_applicable = true,
                                      master_branch = master,
                                      master_branch_is_default_branch = false,
                                      check_license = check_license,
                                      public_registries = public_registries)
                        @info "Running integration test for " test_number master_dir feature_dir public_dir title pass
                        if pass
                            run_thunk()
                        else
                            @test_throws(RegistryCI.AutoMerge.AutoMergeGuidelinesNotMet, run_thunk())
                        end

                    end
                    withenv("AUTOMERGE_GITHUB_TOKEN" => TEST_USER_GITHUB_TOKEN,
                            "TRAVIS_BRANCH" => master,
                            "TRAVIS_BUILD_DIR" => build_dir,
                            "TRAVIS_EVENT_TYPE" => "cron",
                            "TRAVIS_PULL_REQUEST" => "false",
                            "TRAVIS_PULL_REQUEST_SHA" => "",
                            "TRAVIS_REPO_SLUG" => AUTOMERGE_INTEGRATION_TEST_REPO) do
                        sleep(1)
                        AutoMerge.run(;
                                      merge_new_packages = true,
                                      merge_new_versions = true,
                                      new_package_waiting_period = Minute(typemax(Int32)),
                                      new_jll_package_waiting_period = Minute(typemax(Int32)),
                                      new_version_waiting_period = Minute(typemax(Int32)),
                                      new_jll_version_waiting_period = Minute(typemax(Int32)),
                                      registry = AUTOMERGE_INTEGRATION_TEST_REPO,
                                      authorized_authors = String[whoami],
                                      authorized_authors_special_jll_exceptions = String[whoami],
                                      error_exit_if_automerge_not_applicable = true,
                                      master_branch = master,
                                      master_branch_is_default_branch = false)
                        sleep(1)
                        AutoMerge.run(;
                                      merge_new_packages = true,
                                      merge_new_versions = true,
                                      new_package_waiting_period = Minute(0),
                                      new_jll_package_waiting_period = Minute(0),
                                      new_version_waiting_period = Minute(0),
                                      new_jll_version_waiting_period = Minute(0),
                                      registry = AUTOMERGE_INTEGRATION_TEST_REPO,
                                      authorized_authors = String[whoami],
                                      authorized_authors_special_jll_exceptions = String[whoami],
                                      error_exit_if_automerge_not_applicable = true,
                                      master_branch = master,
                                      master_branch_is_default_branch = false)
                    end
                end
            end
        end
    end
end


text:
halve(x::Integer) = x >> one(x)
double(x::Integer) = Int8(2) * x
even(x::Integer) = x & 1 != 1

max_stars_repo_name: ethansaxenian/RosettaDecode

max_issues_repo_name: ethansaxenian/RosettaDecode

max_forks_repo_name: ethansaxenian/RosettaDecode


text:
#= This file is intended to be a drop-in replacement for file rankcorr.jl in the StatsBase
package, except that this file does not contain the code for Spearman's correlation in the
first 116 lines of that file. =#

#######################################
#
#   Kendall correlation
#
#######################################

"""
    corkendall(x, y=x)
Compute Kendall's rank correlation coefficient, τ. `x` and `y` must both be either
matrices or vectors.
"""
function corkendall(x::Union{RealVector,RealOrMissingVector},
                    y::Union{RealVector,RealOrMissingVector};
                    skipmissing::Symbol=:undefined)
    length(x) == length(y) || throw(DimensionMismatch("Vectors must have same length"))
    x, y = handlecompletemissings(x, y, skipmissing)
    ck!(copy(x), copy(y))
end

#= It is idiosyncratic that this method returns a vector, not a matrix, i.e. not consistent
with Statistics.cor or corspearman. But fixing that is a breaking change. =#
function corkendall(X::Union{RealMatrix,RealOrMissingMatrix},
                    y::Union{RealVector,RealOrMissingVector};
                    skipmissing::Symbol=:undefined)
    size(X, 1) == length(y) ||
        throw(DimensionMismatch("X and y have inconsistent dimensions"))
    X, y = handlecompletemissings(X, y, skipmissing)
    n = size(X, 2)
    permy = sortperm(y)
    sortedy = y[permy]
    return([ck_sorted!(copy(sortedy), X[:,i], permy) for i in 1:n])
end

function corkendall(x::Union{RealVector,RealOrMissingVector},
                    Y::Union{RealMatrix,RealOrMissingMatrix};
                    skipmissing::Symbol=:undefined)
    size(Y, 1) == length(x) ||
        throw(DimensionMismatch("x and Y have inconsistent dimensions"))
    x, Y = handlecompletemissings(x, Y, skipmissing)
    n = size(Y, 2)
    permx = sortperm(x)
    sortedx = x[permx]
    return(reshape([ck_sorted!(copy(sortedx), Y[:,i], permx) for i in 1:n], 1, n))
end

function corkendall(X::Union{RealMatrix,RealOrMissingMatrix};
                    skipmissing::Symbol=:undefined)
    X = handlecompletemissings(X, skipmissing)
    n = size(X, 2)
    C = Matrix{Float64}(I, n, n)
    for j = 2:n
        permx = sortperm(X[:,j])
        sortedx = X[:,j][permx]
        for i = 1:j - 1
            C[i,j] = C[j,i] = ck_sorted!(sortedx, X[:,i], permx)
        end
    end
    return C
end

function corkendall(X::Union{RealMatrix,RealOrMissingMatrix},
                    Y::Union{RealMatrix,RealOrMissingMatrix};
                    skipmissing::Symbol=:undefined)
    X, Y = handlecompletemissings(X, Y, skipmissing)
    nr = size(X, 2)
    nc = size(Y, 2)
    C = Matrix{Float64}(undef, nr, nc)
    for j = 1:nr
        permx = sortperm(X[:,j])
        sortedx = X[:,j][permx]
        for i = 1:nc
            C[j,i] = ck_sorted!(sortedx, Y[:,i], permx)
        end
    end
    return C
end

# Knight, William R. “A Computer Method for Calculating Kendall's Tau with Ungrouped Data.”
# Journal of the American Statistical Association, vol. 61, no. 314, 1966, pp. 436–439.
# JSTOR, www.jstor.org/stable/2282833.
"""
    ck_sortedshuffled!(x::RealVector, y::RealVector)
Kendall correlation between two vectors but this function omits the initial sorting of
arguments. So calculating Kendall correlation between `x` and `y` is a three stage process:
a) sort `x` to get `sortedx`; b) apply the same permutation to `y` to get `shuffledy`;
c) call this function on `sortedx` and `shuffledy`.
"""
function ck_sortedshuffled!(x::RealVector, y::RealVector)
    if any(isnan, x) || any(isnan, y) return NaN end
    n = length(x)

    # Use widen to avoid overflows on both 32bit and 64bit
    npairs = div(widen(n) * (n - 1), 2)
    ntiesx = ndoubleties = nswaps = widen(0)
    k = 0

    @inbounds for i = 2:n
        if x[i - 1] == x[i]
            k += 1
        elseif k > 0
            # Sort the corresponding chunk of y, so the rows of hcat(x,y) are
            # sorted first on x, then (where x values are tied) on y. Hence
            # double ties can be counted by calling countties.
            sort!(view(y, (i - k - 1):(i - 1)))
            ntiesx += div(widen(k) * (k + 1), 2) # Must use wide integers here
            ndoubleties += countties(y,  i - k - 1, i - 1)
            k = 0
        end
    end
    if k > 0
        sort!(view(y, (n - k):n))
        ntiesx += div(widen(k) * (k + 1), 2)
        ndoubleties += countties(y, n - k, n)
    end

    nswaps = merge_sort!(y, 1, n)
    ntiesy = countties(y, 1, n)

    # Calls to float below prevent possible overflow errors when
    # length(x) exceeds 77_936 (32 bit) or 5_107_605_667 (64 bit)

    (npairs + ndoubleties - ntiesx - ntiesy - 2 * nswaps) /
     sqrt(float(npairs - ntiesx) * float(npairs - ntiesy))
end

# Auxiliary functions for Kendall's rank correlation
"""
    ck!(x::RealVector, y::RealVector, permx::AbstractVector{<:Integer}=sortperm(x))
Kendall correlation between two vectors `x` and `y`. Third argument `permx` is the
permutation that must be applied to `x` to sort it.
"""
function ck!(x::RealVector, y::RealVector, permx::AbstractVector{<:Integer}=sortperm(x))
    permute!(x, permx)
    permute!(y, permx)
    ck_sortedshuffled!(x, y)
end

function ck!(x::RealOrMissingVector, y::RealOrMissingVector,
            permx::AbstractVector{<:Integer}=sortperm(x))
    permute!(x, permx)
    permute!(y, permx)
    x, y = skipmissingpairs(x, y)
    length(x) >= 2 || return(NaN)
    ck_sortedshuffled!(x, y)
end

"""
    ck_sorted!(sortedx::RealVector, y::RealVector, permx::AbstractVector{<:Integer})
Kendall correlation between two vectors but this function omits the initial sorting of the
first argument. So calculating Kendall correlation between `x` and `y` is a two stage
process: a) sort `x` to get `sortedx`; b) call this function on `sortedx` and `y`, with the
third argument being the permutation that achieved the sorting of `x`.
"""
function ck_sorted!(sortedx::RealVector, y::RealVector, permx::AbstractVector{<:Integer})
    permute!(y, permx)
    ck_sortedshuffled!(sortedx, y)
end
# method for when missings appear, so call skipmissingpairs.
function ck_sorted!(sortedx::RealOrMissingVector, y::RealOrMissingVector,
                    permx::AbstractVector{<:Integer})
    permute!(y, permx)
    sortedx, y = skipmissingpairs(sortedx, y)
    length(sortedx) >= 2 || return(NaN)
    ck_sortedshuffled!(sortedx, y)
end

"""
    countties(x::RealVector, lo::Integer, hi::Integer)
Return the number of ties within `x[lo:hi]`. Assumes `x` is sorted.
"""
function countties(x::AbstractVector, lo::Integer, hi::Integer)
    # Use of widen below prevents possible overflow errors when
    # length(x) exceeds 2^16 (32 bit) or 2^32 (64 bit)
    thistiecount = result = widen(0)
    checkbounds(x, lo:hi)
    @inbounds for i = (lo + 1):hi
        if x[i] == x[i - 1]
            thistiecount += 1
        elseif thistiecount > 0
            result += div(thistiecount * (thistiecount + 1), 2)
            thistiecount = widen(0)
        end
    end

    if thistiecount > 0
        result += div(thistiecount * (thistiecount + 1), 2)
    end
    result
end

# Tests appear to show that a value of 64 is optimal,
# but note that the equivalent constant in base/sort.jl is 20.
const SMALL_THRESHOLD = 64

# merge_sort! copied from Julia Base
# (commit 28330a2fef4d9d149ba0fd3ffa06347b50067647, dated 20 Sep 2020)
"""
    merge_sort!(v::AbstractVector, lo::Integer, hi::Integer, t::AbstractVector=similar(v, 0))
Mutates `v` by sorting elements `x[lo:hi]` using the merge sort algorithm.
This method is a copy-paste-edit of sort! in base/sort.jl, amended to return the bubblesort
distance.
"""
function merge_sort!(v::AbstractVector, lo::Integer, hi::Integer,
                     t::AbstractVector=similar(v, 0))
    # Use of widen below prevents possible overflow errors when
    # length(v) exceeds 2^16 (32 bit) or 2^32 (64 bit)
    nswaps = widen(0)
    @inbounds if lo < hi
        hi - lo <= SMALL_THRESHOLD && return insertion_sort!(v, lo, hi)

        m = midpoint(lo, hi)
        (length(t) < m - lo + 1) && resize!(t, m - lo + 1)

        nswaps = merge_sort!(v, lo,  m, t)
        nswaps += merge_sort!(v, m + 1, hi, t)

        i, j = 1, lo
        while j <= m
            t[i] = v[j]
            i += 1
            j += 1
        end

        i, k = 1, lo
        while k < j <= hi
            if v[j] < t[i]
                v[k] = v[j]
                j += 1
                nswaps += m - lo + 1 - (i - 1)
            else
                v[k] = t[i]
                i += 1
            end
            k += 1
        end
        while k < j
            v[k] = t[i]
            k += 1
            i += 1
        end
    end
    return nswaps
end

# insertion_sort! and midpoint copied from Julia Base
# (commit 28330a2fef4d9d149ba0fd3ffa06347b50067647, dated 20 Sep 2020)
midpoint(lo::T, hi::T) where T <: Integer = lo + ((hi - lo) >>> 0x01)
midpoint(lo::Integer, hi::Integer) = midpoint(promote(lo, hi)...)

"""
    insertion_sort!(v::AbstractVector, lo::Integer, hi::Integer)
Mutates `v` by sorting elements `x[lo:hi]` using the insertion sort algorithm.
This method is a copy-paste-edit of sort! in base/sort.jl, amended to return the bubblesort
distance.
"""
function insertion_sort!(v::AbstractVector, lo::Integer, hi::Integer)
    if lo == hi return widen(0) end
    nswaps = widen(0)
    @inbounds for i = lo + 1:hi
        j = i
        x = v[i]
        while j > lo
            if x < v[j - 1]
                nswaps += 1
                v[j] = v[j - 1]
                j -= 1
                continue
            end
            break
        end
        v[j] = x
    end
    return nswaps
end

"""
    handlecompletemissings(x::AbstractArray,y::AbstractArray,skipmissing::Symbol)
Handles the case of skipmissing = :complete. This is a simpler case than :pairwise, we
merely need to construct new argument(s) for corkendall by calling skipmissingpairs. The
function also validates skipmissing, throwing an error if invalid.
"""
function handlecompletemissings(x::AbstractArray, y::AbstractArray, skipmissing::Symbol)
    if skipmissing == :complete
        if x isa Matrix || y isa Matrix
            return(skipmissingpairs(x, y))
        end
    elseif skipmissing == :pairwise
    elseif skipmissing == :undefined
        if missing isa eltype(x) || missing isa eltype(y)
            throw(ArgumentError("When Missing is an allowed element type"
                * " then keyword argument skipmissing must be either "
                * "`:pairwise` or `:complete`, but got `:$skipmissing`"))
        end
    else
        if missing isa eltype(x) || missing isa eltype(y)
            throw(ArgumentError("keyword argument skipmissing must be either " *
            "`:pairwise` or `:complete`, but got `:$skipmissing`"))
        else
            throw(ArgumentError("keyword argument skipmissing must be either " *
            "`:pairwise`, `:complete` or `:undefined` but got `:$skipmissing`"))
        end
    end
    return(x, y)
end

function handlecompletemissings(x::AbstractArray, skipmissing::Symbol)
    if skipmissing == :complete
        if x isa Matrix
            return(skipmissingpairs(x))
        end
    elseif skipmissing == :pairwise
    elseif skipmissing == :undefined
        if missing isa eltype(x)
            throw(ArgumentError("When Missing is an allowed element type"
                * " then keyword argument skipmissing must be either "
                * "`:pairwise` or `:complete`, but got `:$skipmissing`"))
        end
    else
        if missing isa eltype(x)
            throw(ArgumentError("keyword argument skipmissing must be either " *
            "`:pairwise` or `:complete`, but got `:$skipmissing`"))
        else
            throw(ArgumentError("keyword argument skipmissing must be either " *
            "`:pairwise`, `:complete` or `:undefined` but got `:$skipmissing`"))
        end
    end
    return(x)
end


text:
@connector function Con(;name)
    sts = @variables begin
        P(t)
        Q(t), [connect = Flow]
    end
    ODESystem(Equation[], t, sts, []; name=name)
end

function CompliantCompartment(; name, PCargs...)

    @named pressurizedCompartment = PressurizedCompartment(;PCargs...)

    @named in = Con()
    @named out = Con()

    sts = @variables V(t) P(t) dV(t)

    # Flow is positive into the component.
    eqs = [
        D(V) ~ dV,
        dV ~ in.Q + out.Q,
        in.P ~ P,
        out.P ~ P
    ]

    pressurizedCompliantCompartment = extend(ODESystem(eqs, t, sts, []; name),
                                            pressurizedCompartment)
    compose(pressurizedCompliantCompartment, in, out)
end

"""
This compartment can be thought of as a stiff tube with a fluid that has
mass. Inertia is the resistance to change in velocity of the fluid.

# Arguments
  - `Area`: [m^2] cross-sectional area of the compartment.
  - `V`: [m^3] Volume of the compartment. In the inertia compartment, Volume is a parameter (constant), not a state.
  - `ρ`: [kg/m^3] Density of the fluid.
  - `Q_start`: Initial condition for flow.
"""
function InertiaCompartment(; name,
                               Area::Float64,
                               V::Float64,
                               ρ::Float64 = 1060., # kg/m^3.
                               Q_start::Float64 = 0.)

    @named volumelessComponent = VolumelessComponent()

    @unpack ΔP, Q = volumelessComponent

    @named in = Con()
    @named out = Con()

    ps = @parameters(
        V = V,
        Area = Area,
        ρ = ρ
    )

    # Flow is positive into the component.
    eqs = [
        # Inductance equation.
        # L = ρ * V/Area^2 = rho * length / Area.
        D(Q) ~ ΔP / (ρ * V/Area^2)
    ]

    InertiaCompartment = extend(ODESystem(eqs, t, [], ps; name, defaults = [Q => Q_start]),
                                            volumelessComponent)
    compose(InertiaCompartment, in, out)
end

function VolumelessComponent(; name)
    @named in = Con()
    @named out = Con()

    sts = @variables Q(t) ΔP(t)

    # Flow is positive into the component.
    eqs = [
        0 ~ in.Q + out.Q,
        Q ~ in.Q,
        ΔP ~ in.P - out.P
    ]

    compose(ODESystem(eqs, t, sts, []; name), in, out)
end

"""
# Arguments
  - `ext_pressure` and `tm_pressure`: If set to a number{Float64},
  this will fix the pressure to that value.
  If set to `"free"`, no equation is made for controling
  the state (this should be done in an inheriting system (e.g. `CompliantCompartment()`)).

"""
function PressurizedCompartment(;name, ext_pressure = 0., tm_pressure = "free")

    sts = @variables P(t) P_ext(t) P_tm(t)

    eqs = [
        P ~ P_ext + P_tm
    ]

    if tm_pressure != "free"
        @assert typeof(tm_pressure) == Float64

        push!(eqs, P_tm ~ tm_pressure)
    end

    if ext_pressure != "free"
        @assert typeof(ext_pressure) == Float64

        push!(eqs, P_ext ~ ext_pressure)
    end

    ODESystem(eqs, t, sts, []; name)
end
PressurisedCompartment = PressurizedCompartment