text:
using DataFrames
using LinearAlgebra
using Statistics: std, var, mean, cov
using HighFrequencyCovariance
using Random

brownian_corr_matrix = Hermitian([1.0 0.75 0.5 0.0;
                                  0.0 1.0 0.5 0.25;
                                  0.0 0.0 1.0 0.25;
                                  0.0 0.0 0.0 1.0])
assets = [:BARC, :HSBC, :VODL, :RYAL]
assets1 = [:BARC]
assets2 = [:BARC,  :VODL, :RYAL]
assets3 = [:HSBC,  :VODL, :RYAL]
twister = MersenneTwister(1)

ts1, true_covar, micro_noise, update_rates = generate_random_path(4, 2000; brownian_corr_matrix = brownian_corr_matrix, assets = assets, vols = [0.02,0.03,0.04,0.05], twister = deepcopy(twister), syncronous = true)
ts2, true_covar, micro_noise, update_rates = generate_random_path(4, 2000; brownian_corr_matrix = brownian_corr_matrix, assets = assets, vols = [0.02,0.03,0.04,0.05], twister = deepcopy(twister), syncronous = true)
ts3, true_covar, micro_noise, update_rates = generate_random_path(4, 2000; assets = assets, vols = [0.02,0.03,0.04,0.05], twister = deepcopy(twister), syncronous = true)

iscloser(a,b) = (a.Correlation_error + a.Volatility_error < b.Correlation_error + b.Volatility_error)

# Simple Volatility
simple_vol = simple_volatility(ts1)
all(values(simple_vol) .< 0.1)

# Preav Convergence
preav_estimate1 = preaveraged_covariance(ts1, assets1)
preav_estimate2 = preaveraged_covariance(ts2, assets2)
preav_estimate3 = preaveraged_covariance(ts3, assets3)

# Testing the combination of multiple covariance matrices.
vector_of_covars = [preav_estimate1, preav_estimate2, preav_estimate3]
combined = combine_covariance_matrices(vector_of_covars, [1,2,3], [3,2,1])
isnan(get_correlation(combined, :BARC, :HSBC)) # As they never appear together.

# Testing order does not matter
preav_estimate4 = preaveraged_covariance(ts3, assets)
vector_of_covars = [preav_estimate1, preav_estimate2, preav_estimate4]
combined = combine_covariance_matrices(vector_of_covars, [1,2,3], [3,2,1])
vector_of_covars = [preav_estimate4, preav_estimate2, preav_estimate1]
combined2 = combine_covariance_matrices(vector_of_covars, [3,2,1], [1,2,3])
calculate_mean_abs_distance(combined, combined2).Correlation_error .< 10*eps()
calculate_mean_abs_distance(combined, combined2).Volatility_error .< 10*eps()


text:
isdir("bin") || mkdir("bin")
cd("ratip2013-angular-coefficients")
run(`make`)
##x cd("../grasp2013-nucpot")
##x run(`make`)
cd("..")

ext: jl


text:
# Sequence rules over the [`SimData`](@ref) object,
# calling [`maprule!`](@ref) for each individual `Rule`.
function sequencerules!(simdata::AbstractSimData)
    sequencerules!(simdata, rules(simdata))
end
function sequencerules!(simdata::AbstractSimData, rules::Tuple)
    # Run the first rules
    simdata = maprule!(simdata, rules[1])
    # Run the rest of the rules recursively
    sequencerules!(simdata, tail(rules))
end
sequencerules!(simdata::AbstractSimData, rules::Tuple{}) = simdata


text:
function PlotStems(label_id, x::AbstractArray{<:Real}, y::AbstractArray{<:Real}, args...)
    return PlotStems(label_id, promote(x, y), args...)
end

function PlotStems(values::AbstractArray{T}; count::Integer=length(values),
                   label_id::String="", y_ref=0.0, xscale=1.0, x0=0.0, offset::Integer=0,
                   stride::Integer=1) where {T<:ImPlotData}
    return PlotStems(label_id, values, count, y_ref, xscale, x0, offset, stride * sizeof(T))
end

function PlotStems(values::AbstractArray{T}; count::Integer=length(values),
                   label_id::String="", y_ref=0.0, xscale=1.0, x0=0.0, offset::Integer=0,
                   stride::Integer=1) where {T<:Real}
    return PlotStems(label_id, Float64.(values), count, y_ref, xscale, x0, offset,
                     stride * sizeof(Float64))
end

function PlotStems(x::AbstractArray{T}, y::AbstractArray{T};
                   count::Integer=min(length(x), length(y)), label_id::String="", y_ref=0.0,
                   offset::Integer=0, stride::Integer=1) where {T<:ImPlotData}
    return PlotStems(label_id, x, y, count, y_ref, offset, stride * sizeof(T))
end

function PlotStems(x::AbstractArray{T}, y::AbstractArray{T};
                   count::Integer=min(length(x), length(y)), label_id::String="", y_ref=0.0,
                   offset::Integer=0, stride::Integer=1) where {T<:Real}
    return PlotStems(label_id, Float64.(x), Float64.(y), count, y_ref, offset,
                     stride * sizeof(Float64))
end

function PlotStems(x::AbstractArray{T1}, y::AbstractArray{T2};
                   kwargs...) where {T1<:Real,T2<:Real}
    return PlotBars(promote(x, y)...; kwargs...)
end


text:
# ---
# title: 1160. Find Words That Can Be Formed by Characters
# id: problem1160
# author: Tian Jun
# date: 2020-10-31
# difficulty: Easy
# categories: Array, Hash Table
# link: <https://leetcode.com/problems/find-words-that-can-be-formed-by-characters/description/>
# hidden: true
# ---
#
# You are given an array of strings `words` and a string `chars`.
#
# A string is _good_  if it can be formed by characters from `chars` (each
# character can only be used once).
#
# Return the sum of lengths of all good strings in `words`.
#
#
#
# **Example 1:**
#
#
#
#     Input: words = ["cat","bt","hat","tree"], chars = "atach"
#     Output: 6
#     Explanation:
#     The strings that can be formed are "cat" and "hat" so the answer is 3 + 3 = 6.
#
#
# **Example 2:**
#
#
#
#     Input: words = ["hello","world","leetcode"], chars = "welldonehoneyr"
#     Output: 10
#     Explanation:
#     The strings that can be formed are "hello" and "world" so the answer is 5 + 5 = 10.
#
#
#
#
# **Note:**
#
#   1. `1 <= words.length <= 1000`
#   2. `1 <= words[i].length, chars.length <= 100`
#   3. All strings contain lowercase English letters only.
#
#
## @lc code=start
using LeetCode

## add your code here:
## @lc code=end


text:
using Stheno, Plots, Random, Statistics
using Stheno: @model

# Sample a few times from the prior.
@model model(α, β, l) = GP(eq(α=α, β=β, l=l));
f = model(1, 15, 1.3);
x = range(-3.0, stop=3.0, length=1000);
scatter(x, rand(f(x, 1e-3), 10);
    markersize=0.1,
    markershape=:xcross,
    markerstrokewidth=0.0,
)

# Compute the logpdf of a sample.
σ² = 1e-3;
y = rand(f(x, σ²));
logpdf(f(x, σ²), y)

# Compute + plot the posterior distribution over the process.
f′ = f | (f(x, σ²) ← y);

x_plot = range(-4.0, stop=4.0, length=1000);
f′_xp = f′(x_plot, 1e-6);
μ′_xp, σ′_xp = mean.(marginals(f′_xp)), std.(marginals(f′_xp));

plt = plot(legend=nothing);
scatter!(plt, x, y;
    markersize=0.1,
    markershape=:xcross,
    markerstrokewidth=0.0,
);
plot!(plt, x_plot, mean(f′_xp);
    linewidth=2,
    linecolor=:blue,
);
# plot!(plt, x_plot, mean(f′_xp) .+ 3 .* std.(marginals(f′_xp)));
plot!(plt, x_plot, [μ′_xp μ′_xp];
    linewidth=0.0,
    linecolor=:blue,
    fillrange=[μ′_xp .- 3 .* σ′_xp, μ′_xp .+ 3 * σ′_xp],
    fillalpha=0.3,
    fillcolor=:blue,
);

# Compute posterior distribution over _observations_.
y′_xp = f′(x_plot, σ²);
μ′_xp, σ′_xp = mean.(marginals(y′_xp)), std.(marginals(y′_xp));
plot!(plt, x_plot, [μ′_xp μ′_xp];
    linewidth=0.0,
    linecolor=:blue,
    fillrange=[μ′_xp .- 3 .* σ′_xp, μ′_xp .+ 3 * σ′_xp],
    fillalpha=0.3,
    fillcolor=:red,
);

display(plt);


text:
export getMeasurements, getMeasurementsFD, getMeasurementsLowLevel

function measDataConv(f::MPIFile, args...)
  data = measDataTD(f, args...)

  if eltype(data) <: Integer
    data = map(Float32, data)
  end
  a = rxDataConversionFactor(f)
  if a!=nothing
    for d=1:size(data,2)
      slice = view(data,:,d,:,:)
      rmul!(slice, a[1,d])
      slice .+= a[2,d]
    end
  end
  return data
end

hannWindow(M) = (1.0 .- cos.(2*π/(M-1)*(0:(M-1))))/(M-1)*M

function measDataSpectralLeakageCorrectedSinglePatch(f::MPIFile, frames)
  @debug "Apply Spectral Cleaning"
  numTimePoints = rxNumSamplingPoints(f)
  numReceivers = rxNumChannels(f)
  numFrames = acqNumFrames(f)

  data = zeros(Float32, numTimePoints, numReceivers, 1, length(frames))

  window3 = hannWindow(numTimePoints*3)
  window2 = hannWindow(numTimePoints*2)

  for (i,fr) in enumerate(frames)
    for r in 1:numReceivers
      if fr==1
        tmp = measDataConv(f, fr:fr+1, 1, r)
        data[:,r,1,i] = 1/2 * (tmp[:,1,1,1] .* window2[1:numTimePoints]
                          +  tmp[:,1,1,2] .* window2[1+numTimePoints:2*numTimePoints]
                          );
      elseif fr==numFrames
        tmp = measDataConv(f, fr-1:fr, 1, r)
        data[:,r,1,i] = 1/2 * (tmp[:,1,1,1] .* window2[1:numTimePoints]
                          +    tmp[:,1,1,2] .* window2[1+numTimePoints:2*numTimePoints]
                          );
      else
        tmp = measDataConv(f, fr-1:fr+1, 1, r)
        data[:,r,1,i] = 1/3 * (tmp[:,1,1,1] .* window3[1:numTimePoints]
                          +  tmp[:,1,1,2] .* window3[1+numTimePoints:2*numTimePoints]
                          +  tmp[:,1,1,3] .* window3[1+2*numTimePoints:3*numTimePoints]
                          );
      end
    end
  end
  return data
end

function measDataSpectralLeakageCorrectedMultiPatch(f::MPIFile, frames, periods)
  @debug "Apply Spectral Cleaning"
  numTimePoints = rxNumSamplingPoints(f)
  numReceivers = rxNumChannels(f)
  numFrames = acqNumFrames(f)
  numPeriods = acqNumPeriodsPerFrame(f)

  data = zeros(Float32, numTimePoints, numReceivers, length(periods), length(frames))

  window3 = hannWindow(numTimePoints*3)
  window2 = hannWindow(numTimePoints*2)

  for (i,fr) in enumerate(frames)
    for (p,pe) in enumerate(periods)
      for r in 1:numReceivers
        if pe == 1
          tmp = measDataConv(f, fr, pe:pe+1, r)
          data[:,r,p,i] = 1/2 * (tmp[:,1,1,1] .* window2[1:numTimePoints]
                            +  tmp[:,1,2,1] .* window2[1+numTimePoints:2*numTimePoints]
                            );
        elseif pe == numPeriods
          tmp = measDataConv(f, fr, pe-1:pe, r)
          data[:,r,p,i] = 1/2 * (tmp[:,1,1,1] .* window2[1:numTimePoints]
                            +    tmp[:,1,2,1] .* window2[1+numTimePoints:2*numTimePoints]
                            );
        else
          tmp = measDataConv(f, fr, pe-1:pe+1, r)
          data[:,r,p,i] = 1/3 * (tmp[:,1,1,1] .* window3[1:numTimePoints]
                            +  tmp[:,1,2,1] .* window3[1+numTimePoints:2*numTimePoints]
                            +  tmp[:,1,3,1] .* window3[1+2*numTimePoints:3*numTimePoints]
                            );
        end
      end
    end
  end
  return data
end

function measDataSpectralLeakageCorrected(f::MPIFile, frames, periods)
  if acqNumPeriodsPerFrame(f) == 1
    return measDataSpectralLeakageCorrectedSinglePatch(f, frames)
  else
    return measDataSpectralLeakageCorrectedMultiPatch(f, frames, periods)
  end
end

function measDataLowLevel(f::MPIFile, args...; spectralLeakageCorrection=false)
  if !spectralLeakageCorrection || measIsSpectralLeakageCorrected(f) ||
      acqNumFrames(f) == 1
     tmp = measDataConv(f, args...)
  else
     tmp = measDataSpectralLeakageCorrected(f, args...)
  end
  return tmp
end

function returnasreal(u::AbstractArray{Complex{T}}) where {T}
  return copy(reshape(reinterpret(T,vec(u)),tuple(size(u,1)*2,size(u)[2:end]...)))
end
returnasreal(u::AbstractArray{T}) where {T<:Real} = u

function getAveragedMeasurements(f::MPIFile; frames=1:acqNumFrames(f),
            numAverages=1,  periods=1:acqNumPeriodsPerFrame(f),
            averagePeriodsPerPatch=false, kargs...)

  @debug "frequency and frame selection" rxNumSamplingPoints(f) rxNumChannels(f) acqNumFrames(f)

  if numAverages == 1
    data = measDataLowLevel(f, frames, periods; kargs...)
  else
    nFrames = length(frames)
    nBlocks = ceil(Int, nFrames / numAverages)

    if rem(nFrames, numAverages) != 0
      @warn "numAverages no integer divisor of nFrames.
              Last Block will be averaged over less than $numAverages Frames."
    end

    data = zeros(Float32, rxNumSamplingPoints(f), rxNumChannels(f), acqNumPeriodsPerFrame(f), nBlocks)

    for i = 1:nBlocks
      index1 = 1 + (i-1)*numAverages
      index2 = min( index1 + numAverages-1, nFrames) # ensure that modulo is taken into account

      tmp = measDataLowLevel(f, frames[index1:index2], periods; kargs...)
      data[:,:,:,i] = mean(tmp,dims=4)
    end
  end

  if averagePeriodsPerPatch
    if periods != 1:acqNumPeriodsPerFrame(f)
      error("Option averagePeriodsPerPatch can only be used when all periods are selected")
    end
    data_ = reshape(data, rxNumSamplingPoints(f), rxNumChannels(f),
                          acqNumPeriodsPerPatch(f), acqNumPatches(f), size(data,4))
    dataAv = mean(data,dims=3)

    return reshape(dataAv, rxNumSamplingPoints(f), rxNumChannels(f),
                           acqNumPatches(f), size(data,4))
  else
    return data
  end
end

"""
  getMeasurements(f, [neglectBGFrames]; kargs...) => Array{Float32,4}

Load the measurement data in time domain

Supported keyword arguments:
* frames
* bgCorrection
* interpolateBG
* tfCorrection
* sortFrames
* numAverages
* spectralLeakageCorrection
"""
function getMeasurements(f::MPIFile, neglectBGFrames=true;
      frames=neglectBGFrames ? (1:acqNumFGFrames(f)) : (1:acqNumFrames(f)),
      bgCorrection=false, interpolateBG=false, tfCorrection=rxHasTransferFunction(f),
      sortFrames=false, kargs...)

  if neglectBGFrames
    idx = measFGFrameIdx(f)

    data = getAveragedMeasurements(f; frames=idx[frames], kargs...)

    if bgCorrection
      @debug "Applying bg correction ..."
      idxBG = measBGFrameIdx(f)
      dataBG = getAveragedMeasurements(f; frames=idxBG, kargs...)
      if interpolateBG
        blockLen = measBGFrameBlockLengths(measIsBGFrame(f))
        st = 1
        for j=1:length(blockLen)
          dataBG[:,:,:,st:st+blockLen[j]-1] .=
               mean(dataBG[:,:,:,st:st+blockLen[j]-1], dims=4)
          st += blockLen[j]
        end

        dataBGInterp = interpolate(dataBG,
          (NoInterp(), NoInterp(), NoInterp(), BSpline(Linear()))) #OnCell?

        origIndex = idx[frames]
        M = size(data,4)
        K = size(dataBG,4)
        N = M + K
        for m=1:M
          alpha = (origIndex[m]-1)/(N-1)*(K-1)+1
          for k1=1:size(data,1)
            for k2=1:size(data,2)
              for k3=1:size(data,3)
                data[k1,k2,k3,m] -= dataBGInterp(k1,k2,k3,alpha)
              end
            end
          end
        end
      else
        data[:,:,:,:] .-= mean(dataBG, dims=4)
      end
    end

    if sortFrames
      if calibIsMeanderingGrid(f)
        data[:,:,:,:] = data[:,:,:,meanderingFramePermutation(f)]
      end
    end
  else
    if sortFrames
      permJoint = fullFramePermutation(f)
      data = getAveragedMeasurements(f; frames=permJoint, kargs...)
    else
      data = getAveragedMeasurements(f; frames=frames, kargs...)
    end

    if bgCorrection
      idxBG = measBGFrameIdx(f)
      dataBG = getAveragedMeasurements(f; frames=idxBG, kargs...)

      data[:,:,:,:] .-= mean(dataBG, dims=4)
    end
  end

  if tfCorrection && !measIsTFCorrected(f)
    tf = rxTransferFunction(f)
    inductionFactor = rxInductionFactor(f)

    J = size(data,1)
    dataF = rfft(data, 1)
    dataF ./= tf
    if inductionFactor != nothing
      dataF ./= inductionFactor
    end
    data = irfft(dataF,J,1)
  end

  return data
end

"""
  getMeasurementsFD(f, [neglectBGFrames]; kargs...) => Array{ComplexF32,4}

Load the measurement data in frequency domain

Supported keyword arguments:
* frames
* bgCorrection
* interpolateBG
* tfCorrection
* sortFrames
* numAverages
* spectralLeakageCorrection
* loadasreal
* transposed
* frequencies
"""
function getMeasurementsFD(f::MPIFile, args...;
      loadasreal=false, transposed=false, frequencies=nothing,
      tfCorrection=rxHasTransferFunction(f),  kargs...)

  data = getMeasurements(f, args..., tfCorrection=false; kargs...)

  data = rfft(data, 1)

  if tfCorrection && !measIsTFCorrected(f)
    tf = rxTransferFunction(f)
    data ./= tf
  end

  if frequencies != nothing
    # here we merge frequencies and channels
    data = reshape(data, size(data,1)*size(data,2), size(data,3), size(data,4))
    data = data[frequencies, :, :]
  end

  if transposed
    if frequencies != nothing
      data = permutedims(data, [3,1,2])
    else
      data = permutedims(data, [4,1,2,3])
    end
  end

  if loadasreal
    data = returnasreal(data)
  end

  return data
end

function spectralLeakageCorrectedData(dataIn)
  @debug "Apply Spectral Cleaning"

  numTimePoints = size(dataIn,1)
  numFrames = size(dataIn,2)

  dataOut = zeros(Float32, numTimePoints, numFrames)

  window3 = hannWindow(numTimePoints*3)
  window2 = hannWindow(numTimePoints*2)

  for (i,fr) in enumerate(collect(1:numFrames))
      if fr==1
        dataOut[:,i] = 1/2 * (dataIn[:,fr] .* window2[1:numTimePoints]
                          +  dataIn[:,fr+1] .* window2[1+numTimePoints:2*numTimePoints]
                          );
      elseif fr==numFrames
        dataOut[:,i] = 1/2 * (dataIn[:,fr-1] .* window2[1:numTimePoints]
                          +    dataIn[:,fr] .* window2[1+numTimePoints:2*numTimePoints]
                          );
      else
        dataOut[:,i] = 1/3 * (dataIn[:,fr-1] .* window3[1:numTimePoints]
                          +  dataIn[:,fr] .* window3[1+numTimePoints:2*numTimePoints]
                          +  dataIn[:,fr+1] .* window3[1+2*numTimePoints:3*numTimePoints]
                          );
      end
    end
  return dataOut
end


text:
@testset "Genera" begin
  K, a = MaximalRealSubfield(8, "a")
  Kt, t = PolynomialRing(K, "t")
  L, b = number_field(t^2 - gen(K) * t + 1)
  p = prime_decomposition(maximal_order(K), 2)[1][1]
  G = local_genera_hermitian(L, p, 4, 2, 4, true)
  @test length(G) == 15
  for i in 1:length(G)
    @test rank(G[i]) == 4
    @test representative(G[i]) in G[i]
  end

  p = prime_decomposition(maximal_order(K), 17)[1][1]
  G = local_genera_hermitian(L, p, 5, 5, 5)
  @test length(G) == 7
  for i in 1:length(G)
    @test rank(G[i]) == 5
    @test representative(G[i]) in G[i]
  end

  K, a = MaximalRealSubfield(8, "a")
  Kt, t = K["t"]
  L, b = number_field(t^2 - gen(K) * t + 1)
  p = prime_decomposition(maximal_order(K), 2)[1][1]
  l =  [[(0, 3, 1, 0), (4, 1, 1, 2)],
       [(0, 3, -1, 0), (4, 1, 1, 2)],
       [(0, 3, 1, 0), (4, 1, -1, 2)],
       [(0, 3, -1, 0), (4, 1, -1, 2)],
       [(0, 2, 1, 0), (2, 2, 1, 1)],
       [(0, 2, -1, 0), (2, 2, 1, 1)],
       [(0, 2, 1, 1), (2, 2, 1, 1)],
       [(0, 2, 1, 0), (2, 2, 1, 2)],
       [(0, 2, 1, 1), (2, 2, -1, 1)],
       [(0, 2, -1, 0), (2, 2, 1, 2)],
       [(0, 2, 1, 1), (2, 2, 1, 2)],
       [(0, 1, 1, 0), (1, 2, 1, 1), (2, 1, 1, 1)],
       [(0, 1, -1, 0), (1, 2, 1, 1), (2, 1, 1, 1)],
       [(1, 4, 1, 1)],
       [(1, 4, -1, 1)]]
  Gs = map(x -> genus(HermLat, L, p, x), l)
  myG = local_genera_hermitian(L, p, 4, 2, 4)
  @test length(Gs) == length(myG)
  @test all(x -> x in Gs, myG)
  @test all(x -> x in myG, Gs)

  # Representatives

  K, a = NumberField(x - 1, "a")
  Kt, t = K["t"]
  E, b = NumberField(t^2 + 1, "b")
  p = prime_decomposition(maximal_order(K), 2)[1][1]
  G = genus(HermLat, E, p, [(0, 3, -1, 0)])
  L = representative(G)
  @test length(Hecke.genus_representatives(L)) == 1
end


text:
function vizualize(ds::Observable{Dataset}, p ; linewidth=5, markersize=20)
    f = Figure(resolution=(1500,1200))
    nts, npd = length(ds[].timeseries), length(ds[].ponctual)
    nfigrows = max(nts, npd)
    axs = Matrix{Axis}(undef,(nfigrows,2))
    for (i,tsname) in enumerate(name.(ds[].timeseries))
        name_sym = Symbol(tsname)
        ts = @lift $ds.timeseries[i]
        target = ts[].target
        plot_props = getproperty(p.plot,name_sym)
        axs[i,1] = Axis(f[i, 1]; plot_props.axis...)
        draw_axis(axs[i,1], ts, p ; linewidth, markersize)
        # if (target == [:ϵ̇c]) | (target == [:ϵ̇])
        #     xlims!(axs[i,1], low = 1e-10)
        #     ylims!(axs[i,1], low = 1e-10)
        #     axs[i,1].xscale = log10
        #     axs[i,1].yscale = log10
        # end
    end
    for (i,pdname) in enumerate(name.(ds[].ponctual))
        name_sym = Symbol(pdname)
        plot_props = getproperty(p.plot,name_sym)
        axs[i,2] = Axis(f[i, 2]; plot_props.axis...)
        pd = @lift $ds.ponctual[i]
        draw_axis(axs[i,2], pd, p ; linewidth, markersize)
    end
    assigned_axs = [isassigned(axs,i) for i in eachindex(axs)]
    on(x->autolimits!.(axs[assigned_axs]),ds)
    display(f)
    return f
end

function vizualize(ds::Observable{<:DatasetType}, p ; linewidth=5, markersize=20)
    f = Figure(resolution=(1500,1200))
    name_sym = name(ds[])
    plot_props = getproperty(p.plot,Symbol(name_sym))
    ax = Axis(f[1, 1]; plot_props.axis...)
    draw_axis(ax, ds, p ; linewidth, markersize)
    on(x -> autolimits!(ax), ds)
    display(f)
    f, ax
end

function draw_axis(ax, ds::Observable{TimeseriesDataset}, p ; linewidth=5, markersize=20)
    dg = @lift group($ds)
    name_sym = Symbol(name(ds[]))
    cols_to_plot = getproperty(p.plot,name_sym).syms
    l = s = nothing
    for i in 1:length(dg[])

        x_vec = @lift abs.($dg[i][:, cols_to_plot[1]])
        y_vec = @lift abs.($dg[i][:, cols_to_plot[2]])
        x_sim_vec = cols_to_plot[1] ∈ ds[].target ? @lift(abs.($dg[i][:, string.(cols_to_plot[1])*"_sim"])) : x_vec
        y_sim_vec = cols_to_plot[2] ∈ ds[].target ? @lift(abs.($dg[i][:, string.(cols_to_plot[2])*"_sim"])) : y_vec

        s = scatter!(ax, x_vec, y_vec;
            #marker = markers[j],
            markersize,
            color = Cycled(i)
        )

        l = lines!(ax, x_sim_vec, y_sim_vec;
        color = :black, linewidth
        )

    end
    #xlabel!(ax,x_name)
    #ylabel!(ax, string.(ds[].target)...)
    ax.xlabel = string(cols_to_plot[1])
    ax.ylabel = string(cols_to_plot[2])
    return l, s
end

function draw_axis(ax, ds::Observable{PonctualDataset}, p ; group_var=:pc, linewidth=5, markersize=20)
    data = @lift $ds.data
    #colnames = names(data[][!,Not(ds[].target)])
    name_sym = Symbol(name(ds[]))
    cols_to_plot = getproperty(p.plot,name_sym).syms
    # if length(ds[].target) == 1
    #     x_name = colnames[argmax(length.(unique.(eachcol(data[][!,Not(ds[].target)]))))]
    #     y_name = string(ds[].target[1])
    # else
    #     x_name = string(ds[].target[1])
    #     y_name = string(ds[].target[2])
    # end
    dg = @lift groupby($data, group_var)
    for i in 1:length(dg[])
        g = @lift $dg[i]
        x_vec = @lift abs.($g[:, cols_to_plot[1]])
        y_vec = @lift abs.($g[:, cols_to_plot[2]])
        x_sim_vec = cols_to_plot[1] ∈ ds[].target ? @lift(abs.($g[:, string.(cols_to_plot[1])*"_sim"])) : x_vec
        y_sim_vec = cols_to_plot[2] ∈ ds[].target ? @lift(abs.($g[:, string.(cols_to_plot[2])*"_sim"])) : y_vec

        # x_vec = isnothing(xvar) ? @lift(abs.($g[:,x_name])) : @lift(abs.($g[!,xvar]))
        # if length(ds[].target) == 1
        #     x_vec_sim = isnothing(xvar) ? @lift(abs.($g[:,x_name])) : @lift(abs.($g[!,xvar]))
        # else
        #     x_vec_sim = isnothing(xvar) ? @lift(abs.($g[:,x_name*"_sim"])) : @lift(abs.($g[!,xvar]))
        # end

        # target_sim_vec = @lift abs.($g[:, y_name*"_sim"])
        # target_vec = @lift abs.($g[:, y_name])

        s = scatter!(ax, x_vec, y_vec;
            #marker = markers[j],
            markersize,
            #color = Cycled(i)
        )
        l = scatterlines!(ax, x_sim_vec, y_sim_vec;
        color = :black, linewidth, markersize
        )
    end
    ax.xlabel = string(cols_to_plot[1])
    ax.ylabel = string(cols_to_plot[2])
    #(x_name ∈ ["ϵ̇", "ϵ̇_dev"]) && (ax.xscale = log10)
    return nothing
end

max_stars_repo_name: Leooop/ParametersEstimator

max_issues_repo_name: Leooop/ParametersEstimator

max_forks_repo_name: Leooop/ParametersEstimator


text:
const HELP_INFO = """

This is MonkeyLang.jl $MONKEY_VERSION by $MONKEY_AUTHOR.

Usage:

monkey run <file> [--vm | --jl]
monkey repl [--vm]
"""

function julia_main(; input::IO = stdin, output::IO = stdout)::Cint
    if length(ARGS) == 0 || ARGS[1] ∉ ["run", "repl"]
        println(output, HELP_INFO)
        return 0
    end

    if ARGS[1] == "run"
        if length(ARGS) < 2
            println(output, "Usage: monkey run <file> [--vm | --jl]")
            return 0
        end

        _, ext = splitext(ARGS[2])
        if ext != ".mo"
            println(output, "ERROR: Only .mo files are supported!")
            return 0
        end

        if !isfile(ARGS[2])
            println(output, "ERROR: File not found!")
            return 0
        end

        code = String(read(open(ARGS[2])))
        if length(ARGS) == 3
            if ARGS[3] == "--vm"
                MonkeyLang.run(code; input, output)
            elseif ARGS[3] == "--jl"
                MonkeyLang.Transpilers.JuliaTranspiler.run(code; input, output)
            else
                println(output, "Usage: monkey run <file> [--vm | --jl]")
                return 0
            end
        else
            MonkeyLang.evaluate(code; input, output)
        end
    else
        MonkeyLang.start_repl(; input, output, use_vm = "--vm" ∈ ARGS)
    end

    return 0
end


text:
Pkg.rm("OpenAIGymAPI")
Pkg.clone("https://github.com/JuliaML/OpenAIGymAPI.jl.git")
Pkg.test("OpenAIGymAPI")


text:
@polyvar x y
@variable(model, t)
S = @set x^2 + y^2 == 1 && y - x^2 == 0.5 && x >=0 && y>=0.5
@constraint(model, x + y >= t, domain = S, maxdegree = 2)
@objective(model, Max, t)
optimize!(model)
println("Solution: $(value(t))")
# Returns t ~ 1.3911

max_stars_repo_name: sums-of-squares/tools

max_issues_repo_name: sums-of-squares/tools

max_forks_repo_name: sums-of-squares/tools


text:
abstract type Container{T} <: AbstractArray{T,1} end

Base.size(var::Container) = size(var.data)
Base.eltype(var::Container{T}) where T = T

Base.getindex(var::Container, i::Int) = var.data[i]
Base.getindex(var::Container, I::Vararg{Int,N}) where N = var.data[I...]
Base.getindex(var::Container, ::Colon) = var.data[:]
Base.getindex(var::Container, kr::AbstractRange) = var.data[kr]

Base.setindex!(var::Container, v, i::Int) = (var.data[i] = v)
Base.setindex!(var::Container, v, I::Vararg{Int,N}) where N = (var.data[I...] = v)
Base.setindex!(var::Container, v, ::Colon) = (var.data[:] .= v)
Base.setindex!(var::Container, v, kr::AbstractRange) = (var.data[kr] .= v)

"""
    struct HFcContainer{T,D,TH,TF,Tc} <: AbstractBuffer{T}
        data::Vector{T}
        H::TH
        F::TF
        c::Tc
    end

A buffer containing data for in-place computations of H,F,c terms.
"""
struct HFcContainer{T,D,TH,TF,Tc} <: Container{T}
    data::Vector{T}
    H::TH
    F::TF
    c::Tc

    function HFcContainer(::Type{T}, D::Integer) where T
        data = zeros(T, D*D+D+1)
        H = reshape( view(data, 1:(D*D)), (D,D) )
        F = view(data, (D*D+1):(D*D+D))
        c = view(data, (D*D+D+1):(D*D+D+1))
        new{T,D,typeof(H),typeof(F),typeof(c)}(data, H, F, c)
    end

    function HFcContainer(data::AbstractArray{T,1}, D) where T
        H = reshape( view(data, 1:(D*D)), (D,D) )
        F = view(data, (D*D+1):(D*D+D))
        c = view(data, (D*D+D+1):(D*D+D+1))
        new{T,D,typeof(H),typeof(F),typeof(c)}(data, H, F, c)
    end
end

function Base.similar(c::HFcContainer{T,D}) where {T,D}
    HFcContainer(similar(c.data), D)
end

function Base.similar(c::HFcContainer{T,D}, ::Type{K}) where {T,D,K}
    HFcContainer(similar(c.data, K), D)
end


text:
using Test

@testset "trained_weights" begin
    include("trained_weights/MNIST.n1.jl")
    include("trained_weights/MNIST.WK17a_linf0.1_authors.jl")
end


text:
#=
Statistics for Unicode character/string types

Copyright 2017 Gandalf Software, Inc., Scott P. Jones
Licensed under MIT License, see LICENSE.md
=#
struct LineCounts
    lines::Int
    chars::Int
    bytes::Int

    LineCounts(lines::Vector) = new(length(lines), sum(length, lines), sum(sizeof, lines))
    LineCounts(l::Integer, c::Integer) = new(Int(l), Int(c), 0)
    LineCounts(l::Integer, c::Integer, b::Integer) = new(Int(l), Int(c), Int(b))
end

LineCounts(old::LineCounts, lines::Vector) = LineCounts(old.lines, old.chars, sum(sizeof, lines))

struct CharTypes
    ascii::Int
    latin::Int
    utf2byte::Int
    ucs2::Int
    utf32::Int
    surr::Int
    invalid::Int
end

struct CharStat
    num::Int
    len::Int
    empty::Int
    minlen::Int
    maxlen::Int
    maxtyp::Int
    chars::CharTypes
    lines::CharTypes
end

maxbit(flags) = sizeof(flags)<<3 - leading_zeros(flags)

function calcstats(lines)
    numlines = length(lines)
    lenlines = sum(length, lines)
    charstat = [0,0,0,0,0,0,0]
    linestat = [0,0,0,0,0,0,0]
    emptylines = 0
    minlen = typemax(Int)
    maxlen = 0
    maxtyp = 0
    for l in lines
        len = length(l)
        minlen = min(minlen, len)
        maxlen = max(maxlen, len)
        flags = 0
        for chr in l
            ch = chr%UInt32
            t = ch <= 0x7f ? 1 :
                ch <= 0xff ? 2 :
                ch <= 0x7ff ? 3 :
                0xd800 <= ch <= 0xdfff ? 6 :
                ch <= 0xffff ? 4 :
                ch <= 0x10ffff ? 5 : 7
            charstat[t] += 1
            flags |= 1<<(t-1)
        end
        if flags == 0
            emptylines += 1
        else
            linestat[maxbit(flags & 0x1f)] += 1
            (flags & (1<<5)) != 0 && (linestat[6] += 1)
            (flags & (1<<6)) != 0 && (linestat[7] += 1)
        end
        maxtyp |= flags
    end
    CharStat(numlines, lenlines, emptylines, minlen, maxlen, maxbit(maxtyp & 0x1f),
             CharTypes(charstat...), CharTypes(linestat...))
end


text:
import Multroot: convolute

conv(a, b) = convolute(a, b)
conv(a::Poly, b::Poly) = conv(p.a, b.a)
function bt01()
#
# Brugnano and Trigiante
#
    p = poly([ones(6);-1;-1;im;im;im;-im;-im;-im;2])
    z = [1 6; -1  2; im  3; -im  3; 2  1];
    p, PolyZeros(z)
 end
function bt02()
#
# Brugnano and Trigiante
#
    p = poly([ones(10);-1;-1;im;-im;2])
    z = [1.0 10; -1 2; im 1; -im 1; 2 1];

    p, PolyZeros(z)
end
function bt03()
#
# Brugnano and Trigiante
#
    p = poly([im*ones(5);-im*ones(5);0.5im*ones(4);-0.5im*ones(4);0.75im;-0.75im])
    z = [im 5; -im  5; 0.5im  4; -0.5im  4; .75im  1; -0.75im  1];

    p, PolyZeros(z)
end
function bt04()
#
# Brugnano and Trigiante
#
    p = poly([1;1;1; -1*[1;1;1;1]; (.5+im)*[1;1;1]; (.5-im)*[1;1;1];
                      .5*(1+im)*[1;1]; .5*(1-im)*[1;1]])

    z = [1 3; -1 4; .5+im 3; .5-im 3; .5*(1+im) 2; .5*(1-im) 2]

    p, PolyZeros(z)
end
function farloi01()
#
# Farmer-Loizou
#
    p = conv([1.0;1;2],[1;1;3]);
    p = conv(p,p); p = conv(p,p);
    z = [-0.50000000000000 + 1.65831239517770im
        -0.50000000000000 - 1.65831239517770im
        -0.50000000000000 + 1.32287565553230im
        -0.50000000000000 - 1.32287565553230im];
    z = [z 4*ones(4)];

    p, PolyZeros(z)
end
function fib(n)
#
#  test polynomial suggested by Goedecker
#
    [-1.0; ones(n)]
end
function fib05()
#
#  test polynomial suggested by Goedecker
#
    n = 5;
    p = fib(n)
    z = [-.6783507129699967-.4585361872731445*im   1;
         -.6783507129699967+.4585361872731445*im   1;
         .19537659464725405-.8488536405462456*im   1;
         .19537659464725405+.8488536405462456*im   1;
         1.9659482366454853                       1];

    p, PolyZeros(z)
end
function fib10()
#
#  test polynomial suggested by Goedecker
#
    n = 10
    p = fib(n)
    z = [-.8990310941957183    1;
         -.7399333995096503-.5167965365015562*im  1;
         -.7399333995096503+.5167965365015562*im  1;
         -.3130259587265545-.8583547482505247*im  1;
         -.3130259587265545+.8583547482505247*im  1;
         .2462366020696891-.9013280059509374*im  1;
         .2462366020696891+.9013280059509374*im  1;
         .7567289869093243-.6038877405652167*im  1;
         .7567289869093243+.6038877405652167*im  1;
         1.999018632710101  1; ]
    p, PolyZeros(z)
end
function fib100()
#
#  test polynomial suggested by Goedecker
#
    n = 100;
    p = fib(n)
    z = [-.9891099736226092;
         -.9871753104718534-.06190245872746059*im;  -.9871753104718534+.06190245872746059*im;
         -.9813788056446001-.1235640063982463*im;  -.9813788056446001+.1235640063982463*im;
         -.9717428840956687-.1847446602565987*im; -.9717428840956687+.1847446602565987*im;
         -.9583048244621342-.2452062905225305*im; -.9583048244621342+.2452062905225305*im;
         -.9411166150274026-.3047135378280153*im; -.9411166150274026+.3047135378280153*im;
         -.9202447528605484-.3630347198287647*im; -.9202447528605484+.3630347198287647*im;
         -.8957699869173070-.4199427234459649*im; -.8957699869173070+.4199427234459649*im;
         -.8677870061092108-.4752158792455410*im; -.8677870061092108+.4752158792455410*im;
         -.8364040735640849-.5286388145285358*im; -.8364040735640849+.5286388145285358*im;
         -.8017426085137990-.5800032817847029*im; -.8017426085137990+.5800032817847029*im;
         -.7639367174532022-.6291089592520593*im; -.7639367174532022+.6291089592520593*im;
         -.7231326764169818-.6757642204274601*im; -.7231326764169818+.6757642204274601*im;
         -.6794883664183336-.7197868694867672*im; -.6794883664183336+.7197868694867672*im;
         -.6331726642844369-.7610048396973263*im; -.6331726642844369+.7610048396973263*im;
         -.5843647913085710-.7992568520396451*im; -.5843647913085710+.7992568520396451*im;
         -.5332536223172150-.8343930313987525*im; -.5332536223172150+.8343930313987525*im;
         -.4800369579227939-.8662754778380540*im; -.4800369579227939+.8662754778380540*im;
         -.4249207628992477-.8947787906289414*im; -.4249207628992477+.8947787906289414*im;
         -.3681183737790589-.9197905428773492*im; -.3681183737790589+.9197905428773492*im;
         -.3098496789278928-.9412117047633791*im; -.3098496789278928+.9412117047633791*im;
         -.2503402745082667-.9589570135917125*im; -.2503402745082667+.9589570135917125*im;
         -.1898205998989927-.9729552890388070*im; -.1898205998989927+.9729552890388070*im;
         -.1285250562957008-.9831496921783526*im; -.1285250562957008+.9831496921783526*im;
         -.06669111238377778-.9894979270705116*im; -.06669111238377778+.9894979270705116*im;
         -.004558401154059704-.9919723839157337*im; -.004558401154059704+.9919723839157337*im;
         .05763218786929409-.9905602230050737*im; .05763218786929409+.9905602230050737*im;
         .1196394164912155-.9852633989536015*im; .1196394164912155+.9852633989536015*im;
         .1812226010328281-.9760986249939283*im; .1812226010328281+.9760986249939283*im;
         .2421425057566853-.9630972774521361*im; .2421425057566853+.9630972774521361*im;
         .3021622210555556-.9463052409577104*im; .3021622210555556+.9463052409577104*im;
         .3610480204611693-.9257826954965046*im; .3610480204611693+.9257826954965046*im;
         .4185701899698170-.9016038471668328*im; .4185701899698170+.9016038471668328*im;
         .4745038225339292-.8738566055392182*im; .4745038225339292+.8738566055392182*im;
         .5286295698615102-.8426422119868863*im; .5286295698615102+.8426422119868863*im;
         .5807343429795305-.8080748254361926*im; .5807343429795305+.8080748254361926*im;
         .6306119525307706-.7702810749337576*im; .6306119525307706+.7702810749337576*im;
         .6780636798327496-.7293995925445015*im; .6780636798327496+.7293995925445015*im;
         .7228987709578736-.6855805456977909*im; .7228987709578736+.6855805456977909*im;
         .7649348495515867-.6389851953915993*im; .7649348495515867+.6389851953915993*im;
         .8039982514449404-.5897855154613771*im; .8039982514449404+.5897855154613771*im;
         .8399242976881979-.5381639173058062*im; .8399242976881979+.5381639173058062*im;
         .8725575452576575-.4843131310830981*im; .8725575452576575+.4843131310830981*im;
         .9017520886609776-.4284362924194566*im; .9017520886609776+.4284362924194566*im;
         .9273720302113955-.3707472629394290*im; .9273720302113955+.3707472629394290*im;
         .9492922833419482-.3114711598054638*im; .9492922833419482+.3114711598054638*im;
         .9673998997347004-.2508449722447981*im; .9673998997347004+.2508449722447981*im;
         .9815960786437149-.1891180053506043*im; .9815960786437149+.1891180053506043*im;
         .9917988816037989-.1265517518780266*im; .9917988816037989+.1265517518780266*im;
         .9979464229745987-.06341873580445710*im; .9979464229745987+.06341873580445710*im;
         2.000000000000000];
    z = [z ones(n)];
    p, PolyZeros(z)
end
function fib15()
#
#  test polynomial suggested by Goedecker
#
    n = 15;
    p = fib(n)
    z = [-.9115849521765970-.1896359497733640*im;
         -.9115849521765970+.1896359497733640*im;
         -.7620785395636261-.5388639564255547*im;
         -.7620785395636261+.5388639564255547*im;
         -.4858038058618722-.8023902123527991*im;
         -.4858038058618722+.8023902123527991*im;
         -.1248130033927766-.9371029958793902*im;
         -.1248130033927766+.9371029958793902*im;
         .2657642744875577-.9183694775911604*im;
         .2657642744875577+.9183694775911604*im;
         .6253568882929779-.7424114494678755*im;
         .6253568882929779+.7424114494678755*im;
         .8931744004970845-.4240640341420445*im;
         .8931744004970845+.4240640341420445*im;
         1.999969475434503];

    z = [z ones(n)];
    p, PolyZeros(z)
end
function fib150()
#
#  test polynomial suggested by Goedecker
#
    n = 150;
    p = fib(n)
    z = [-.9927187559253492;
         -.9918530882118113-.4147913261323188e-1*im;
         -.9918530882118113+.4147913261323188e-1*im;
         -.9892575836821723-.8288617318453414e-1*im;
         -.9892575836821723+.8288617318453414e-1*im;
         -.9849367355771876-.1241491541448754*im;
         -.9849367355771876+.1241491541448754*im;
         -.9788980239954073-.1651963566541725*im;
         -.9788980239954073+.1651963566541725*im;
         -.9711519029549269-.2059564344240104*im;
         -.9711519029549269+.2059564344240104*im;
         -.9617117823114389-.2463585368912788*im;
         -.9617117823114389+.2463585368912788*im;
         -.9505940045641273-.2863324315280673*im;
         -.9505940045641273+.2863324315280673*im;
         -.9378178165898683-.3258086250746164*im;
         -.9378178165898683+.3258086250746164*im;
         -.9234053363550660-.3647184834839322*im;
         -.9234053363550660+.3647184834839322*im;
         -.9073815146632301-.4029943503688435*im;
         -.9073815146632301+.4029943503688435*im;
         -.8897740920050947-.4405696637447960*im;
         -.8897740920050947+.4405696637447960*im;
         -.8706135505866639-.4773790708645437*im;
         -.8706135505866639+.4773790708645437*im;
         -.8499330616190360-.5133585409440989*im;
         -.8499330616190360+.5133585409440989*im;
         -.8277684279621956-.5484454755828406*im;
         -.8277684279621956+.5484454755828406*im;
         -.8041580222231623-.5825788166845422*im;
         -.8041580222231623+.5825788166845422*im;
         -.7791427204169229-.6156991516902623*im;
         -.7791427204169229+.6156991516902623*im;
         -.7527658313064592-.6477488159385299*im;
         -.7527658313064592+.6477488159385299*im;
         -.7250730215458899-.6786719919730538*im;
         -.7250730215458899+.6786719919730538*im;
         -.6961122367582722-.7084148056232615*im;
         -.6961122367582722+.7084148056232615*im;
         -.6659336186869521-.7369254186883471*im;
         -.6659336186869521+.7369254186883471*im;
         -.6345894185664953-.7641541180611353*im;
         -.6345894185664953+.7641541180611353*im;
         -.6021339068661845-.7900534011339715*im;
         -.6021339068661845+.7900534011339715*im;
         -.5686232795658151-.8145780573349891*im;
         -.5686232795658151+.8145780573349891*im;
         -.5341155611300741-.8376852456494883*im;
         -.5341155611300741+.8376852456494883*im;
         -.4986705043541354-.8593345679877696*im;
         -.4986705043541354+.8593345679877696*im;
         -.4623494872592649-.8794881382675888*im;
         -.4623494872592649+.8794881382675888*im;
         -.4252154072232028-.8981106470864277*im;
         -.4252154072232028+.8981106470864277*im;
         -.3873325725358873-.9151694218659986*im;
         -.3873325725358873+.9151694218659986*im;
         -.3487665915767304-.9306344823588107*im;
         -.3487665915767304+.9306344823588107*im;
         -.3095842598151608-.9444785914142191*im;
         -.3095842598151608+.9444785914142191*im;
         -.2698534448415497-.9566773009091463*im;
         -.2698534448415497+.9566773009091463*im;
         -.2296429696409571-.9672089927566195*im;
         -.2296429696409571+.9672089927566195*im;
         -.1890224943274280-.9760549149134032*im;
         -.1890224943274280+.9760549149134032*im;
         -.1480623965618780-.9831992123163535*im;
         -.1480623965618780+.9831992123163535*im;
         -.1068336508819937-.9886289526856903*im;
         -.1068336508819937+.9886289526856903*im;
         -.6540770717811572e-1-.9923341471422300*im;
         -.6540770717811572e-1+.9923341471422300*im;
         -.2385636855484892e-1-.9943077655947946*im;
         -.2385636855484892e-1+.9943077655947946*im;
         .1774833117573815e-1-.9945457468636005*im;
         .1774833117573815e-1+.9945457468636005*im;
         .5933425011689494e-1-.9930470035155510*im;
         .5933425011689494e-1+.9930470035155510*im;
         .1008292607663695-.9898134213981652*im;
         .1008292607663695+.9898134213981652*im;
         .1421613720845280-.9848498538705906*im;
         .1421613720845280+.9848498538705906*im;
         .1832588508751938-.9781641107430547*im;
         .1832588508751938+.9781641107430547*im;
         .2240503421931883-.9697669419505880*im;
         .2240503421931883+.9697669419505880*im;
         .2644649884815679-.9596720160034037*im;
         .2644649884815679+.9596720160034037*im;
         .3044325471292468-.9478958932755882*im;
         .3044325471292468+.9478958932755882*im;
         .3438835061260140-.9344579942165983*im;
         .3438835061260140+.9344579942165983*im;
         .3827491974769262-.9193805625975399*im;
         .3827491974769262+.9193805625975399*im;
         .4209619080219071-.9026886239377113*im;
         .4209619080219071+.9026886239377113*im;
         .4584549872896030-.8844099392981356*im;
         .4584549872896030+.8844099392981356*im;
         .4951629519980681-.8645749546799259*im;
         .4951629519980681+.8645749546799259*im;
         .5310215868002686-.8432167463288874*im;
         .5310215868002686+.8432167463288874*im;
         .5659680408622439-.8203709623267845*im;
         .5659680408622439+.8203709623267845*im;
         .5999409198599586-.7960757609476053*im;
         .5999409198599586+.7960757609476053*im;
         .6328803729932399-.7703717463775087*im;
         .6328803729932399+.7703717463775087*im;
         .6647281746501417-.7433019025433064*im;
         .6647281746501417+.7433019025433064*im;
         .6954278004244162-.7149115259686257*im;
         .6954278004244162+.7149115259686257*im;
         .7249244973086152-.6852481587793065*im;
         .7249244973086152+.6852481587793065*im;
         .7531653480767721-.6543615232056414*im;
         .7531653480767721+.6543615232056414*im;
         .7800993301600754-.6223034591666788*im;
         .7800993301600754+.6223034591666788*im;
         .8056773697373049-.5891278667468168*im;
         .8056773697373049+.5891278667468168*im;
         .8298523923415625-.5548906555452699*im;
         .8298523923415625+.5548906555452699*im;
         .8525793720535097-.5196497029280600*im;
         .8525793720535097+.5196497029280600*im;
         .8738153823196024-.4834648230425889*im;
         .8738153823196024+.4834648230425889*im;
         .8935196525762829-.4463977479358605*im;
         .8935196525762829+.4463977479358605*im;
         .9116536360890624-.4085121210914591*im;
         .9116536360890624+.4085121210914591*im;
         .9281810955442824-.3698735020064003*im;
         .9281810955442824+.3698735020064003*im;
         .9430682136510916-.3305493779521499*im;
         .9430682136510916+.3305493779521499*im;
         .9562837358786358-.2906091758178016*im;
         .9562837358786358+.2906091758178016*im;
         .9677991509344817-.2501242631678211*im;
         .9677991509344817+.2501242631678211*im;
         .9775889111839516-.2091679239550551*im;
         .9775889111839516+.2091679239550551*im;
         .9856306896705420-.1678152916898196*im;
         .9856306896705420+.1678152916898196*im;
         .9919056630215224-.1261432225318404*im;
         .9919056630215224+.1261432225318404*im;
         .9963988013915535-.8423009394365190e-1*im;
         .9963988013915535+.8423009394365190e-1*im;
         .9990991395939178-.4215552183469491e-1*im;
         .9990991395939178+.4215552183469491e-1*im;
         2.000000000000000];

    z = [z ones(n)];
    p, PolyZeros(z)
end
function fib20()
#
#  test polynomial suggested by Goedecker
#
    n = 20;
    p = fib(n)
    z = [-.9473885474109854;
         -.9029230998844866-.2884421994216345*im; -.9029230998844866+.2884421994216345*im;
         -.7734751441751908-.5504669225440185*im; -.7734751441751908+.5504669225440185*im;
         -.5705411148090332-.7619498890625472*im; -.5705411148090332+.7619498890625472*im;
         -.3121525651549263-.9031354642584495*im; -.3121525651549263+.9031354642584495*im;
         -.02129876087300105-.9602891496489760*im; -.02129876087300105+.9602891496489760*im;
         .2760438985065760-.9267417870492216*im; .2760438985065760+.9267417870492216*im;
         .5529956567245064-.8031600819414702*im; .5529956567245064+.8031600819414702*im;
         .7834063973613306-.5969423289081712*im; .7834063973613306+.5969423289081712*im;
         .9416394828514235-.3213776783384489*im; .9416394828514235+.3213776783384489*im;
         1.999999046316589];
    z = [z ones(n)]
    p, PolyZeros(z)
end
function fib30()
#
#  test polynomial suggested by Goedecker
#
    n = 30;
    p = fib(30)
    z = [-.9644254204388659;
         -.9439446810207413-.1984294062681023*im; -.9439446810207413+.1984294062681023*im;
         -.8833406981383851-.3885717512211873*im; -.8833406981383851+.3885717512211873*im;
         -.7850939960122262-.5624744046156477*im; -.7850939960122262+.5624744046156477*im;
         -.6532262984053521-.7128393846266873*im; -.6532262984053521+.7128393846266873*im;
         -.4931368697016981-.8333156264808382*im; -.4931368697016981+.8333156264808382*im;
         -.3113834310599106-.9187504784857605*im; -.3113834310599106+.9187504784857605*im;
         -.1154175582352257-.9653888511274815*im; -.1154175582352257+.9653888511274815*im;
         .8671260421202705e-1-.9710099019085114*im; .8671260421202705e-1+.9710099019085114*im;
         .2866757308842761-.9349927196966469*im; .2866757308842761+.9349927196966469*im;
         .4761639390711853-.8583045796885356*im; .4761639390711853+.8583045796885356*im;
         .6471568601471165-.7434108410885516*im; .6471568601471165+.7434108410885516*im;
         .7920805342279742-.5941279591518842*im; .7920805342279742+.5941279591518842*im;
         .9037637390940364-.4155315939426088*im; .9037637390940364+.4155315939426088*im;
         .9752028356220179-.2142999451221434*im; .9752028356220179+.2142999451221434*im;
         1.999999999068677];

    z = [z ones(n)]
    p, PolyZeros(z)
end
function fib50()
#
#  test polynomial suggested by Goedecker
#
    n = 50;
    p = fib(n)
    z = [-.9784087187483260;
         -.9708270453807137-.1218355210630181*im; -.9708270453807137+.1218355210630181*im;
         -.9481969430968524-.2418020920979234*im; -.9481969430968524+.2418020920979234*im;
         -.9108614275411458-.3580588620158148*im; -.9108614275411458+.3580588620158148*im;
         -.8593864256832723-.4688207429141271*im; -.8593864256832723+.4688207429141271*im;
         -.7945522211328693-.5723852109383657*im; -.7945522211328693+.5723852109383657*im;
         -.7173416644328307-.6671578269195402*im; -.7173416644328307+.6671578269195402*im;
         -.6289253327402240-.7516760773713197*im; -.6289253327402240+.7516760773713197*im;
         -.5306438727753174-.8246311589916703*im; -.5306438727753174+.8246311589916703*im;
         -.4239878085206529-.8848873569443235*im; -.4239878085206529+.8848873569443235*im;
         -.3105751412168900-.9314986982161571*im; -.3105751412168900+.9314986982161571*im;
         -.1921271147696794-.9637225955421171*im; -.1921271147696794+.9637225955421171*im;
         -.07044256696033078-.9810302342100871*im; -.07044256696033078+.9810302342100871*im;
         .5262866000569139e-1-.9831134936221445*im; .5262866000569139e-1+.9831134936221445*im;
         .1752127097133740-.9698882397909325*im; .1752127097133740+.9698882397909325*im;
         .2954384705313263-.9414938807153827*im; .2954384705313263+.9414938807153827*im;
         .4114633335560060-.8982891618708299*im; .4114633335560060+.8982891618708299*im;
         .5214969674899803-.8408443381786473*im; .5214969674899803+.8408443381786473*im;
         .6238218132212410-.7699301944269942*im; .6238218132212410+.7699301944269942*im;
         .7168085997048416-.6865051209865677*im; .7168085997048416+.6865051209865677*im;
         .7989249239308063-.5917030159705396*im; .7989249239308063+.5917030159705396*im;
         .8687357356619744-.4868277941747678*im; .8687357356619744+.4868277941747678*im;
         .9248989976696934-.3733647458567231*im; .9248989976696934+.3733647458567231*im;
         .9661725191057798-.2530209270696782*im; .9661725191057798+.2530209270696782*im;
         .9914691930342274-.1277916638512058*im; .9914691930342274+.1277916638512058*im;
         1.999999999999999];

    z = [z ones(n)]
    p, PolyZeros(z)
end
function fibsq04()
#
#  test polynomial suggested by Goedecker
#
    n = 4;
    p = fib(n)
    p = conv(p,p)
    z = [-.7748041132154339;
         -.07637893113374573-.8147036471703865*im; -.07637893113374573+.8147036471703865*im;
         1.927561975482925]
    z = [z  2*ones(n)]
    p, PolyZeros(z)
end
function fibsq08()
#
#  test polynomial suggested by Goedecker
#
    n = 8;
    p = fib(n)
    p = conv(p,p)
    z = [-.8762862300182460;
         -.6416053894725764-.6063952206057722*im; -.6416053894725764+.6063952206057722*im;
         -.04694032461065018-.9030234661229178*im; -.04694032461065018+.9030234661229178*im;
         .6286732392246423-.7084725692273577*im; .6286732392246423+.7084725692273577*im;
         1.996031179735415];
    z = [z 2*ones(n)]
    p, PolyZeros(z)
end
function fibsq16()
#
#  test polynomial suggested by Goedecker
#  square of Fibocacci polynomial
#
    n = 16;
    p = [-1; ones(n)];
    p = conv(p,p);
    z = [-.9349216524104182;
         -.8673161854070643-.3515310211756676*im;
         -.8673161854070643+.3515310211756676*im;
         -.6736211580446324-.6537256870201036*im;
         -.6736211580446324+.6537256870201036*im;
         -.3799784643132753-.8637046142557905*im;
         -.3799784643132753+.8637046142557905*im;
         -.2606735312578153e-1-.9505412093420680*im;
         -.2606735312578153e-1+.9505412093420680*im;
         .3400699024566086-.8988910323302336*im;
         .3400699024566086+.8988910323302336*im;
         .6677317831101087-.7098987824372997*im;
         .6677317831101087+.7098987824372997*im;
         .9066499318552733-.3989926258132711*im;
         .9066499318552733+.3989926258132711*im;
         1.999984739347944];
    z = [z 2*ones(n)]
    p, PolyZeros(z)
end
function fibsq24()
#
#  test polynomial suggested by Goedecker
#  square of Fibocacci polynomial
#
    n = 24;
    p = fib(n)
    p = conv(p, p)
    z = [-.9558467190743565;
         -.9244083656460181-.2442575811931259*im;
         -.9244083656460181+.2442575811931259*im;
         -.8320677914606713-.4727784060949767*im;
         -.8320677914606713+.4727784060949767*im;
         -.6846250033032301-.6707962641470166*im;
         -.6846250033032301+.6707962641470166*im;
         -.4913433576989534-.8254218906657225*im;
         -.4913433576989534+.8254218906657225*im;
         -.2643730603203030-.9264240880183555*im;
         -.2643730603203030+.9264240880183555*im;
         -.1800084671534713e-1-.9668298465915332*im;
         -.1800084671534713e-1+.9668298465915332*im;
         .2322158170552186-.9432943625623267*im;
         .2322158170552186+.9432943625623267*im;
         .4703449187206489-.8561989307649443*im;
         .4703449187206489+.8561989307649443*im;
         .6808619003199168-.7094507040001145*im;
         .6808619003199168+.7094507040001145*im;
         .8490194865023282-.5100588510854889*im;
         .8490194865023282+.5100588510854889*im;
         .9602996918859324-.2681678908797330*im;
         .9602996918859324+.2681678908797330*im;
         1.999999940395313];
    z = [z 2*ones(n)]
    p, PolyZeros(z)
end
function fibsq32()
#
#  test polynomial suggested by Goedecker
#  square of Fibocacci polynomial
#
    n = 32;
    p = fib(n)
    p = conv(p,p)
    z = [-.9665892387165332;
         -.9485092065867169-.1867222182129376*im;
         -.9485092065867169+.1867222182129376*im;
         -.8949223987510440-.3665688458356083*im;
         -.8949223987510440+.3665688458356083*im;
         -.8077651569138829-.5329094943253877*im;
         -.8077651569138829+.5329094943253877*im;
         -.6901871862102008-.6795950061744879*im;
         -.6901871862102008+.6795950061744879*im;
         -.5464383280299697-.8011754158225671*im;
         -.5464383280299697+.8011754158225671*im;
         -.3817161237254617-.8930914814237271*im;
         -.3817161237254617+.8930914814237271*im;
         -.2019801534661303-.9518321713734013*im;
         -.2019801534661303+.9518321713734013*im;
         -.1374081948052525e-1-.9750513807247308*im;
         -.1374081948052525e-1+.9750513807247308*im;
         .1761678077924798-.9616381291177512*im;
         .1761678077924798+.9616381291177512*im;
         .3608186423377573-.9117356125589683*im;
         .3608186423377573+.9117356125589683*im;
         .5334081767292508-.8267064135873833*im;
         .5334081767292508+.8267064135873833*im;
         .6874310990889808-.7090471025782057*im;
         .6874310990889808+.7090471025782057*im;
         .8167660514846921-.5622778345619681*im;
         .8167660514846921+.5622778345619681*im;
         .9156129648696007-.3909102172601782*im;
         .9156129648696007+.3909102172601782*im;
         .9783492503358520-.2007808902007197*im;
         .9783492503358520+.2007808902007197*im;
         1.999999999767169];
    z = [z 2*ones(n)]
    p, PolyZeros(z)
end
function fibsq48()
#
#  test polynomial suggested by Goedecker
#  square of Fibocacci polynomial
#
    n = 48;
    p = fib(n)
    p = conv(p,p)
    z = [-.9775253106342052;
         -.9693128478789493-.1267362311985826*im;
         -.9693128478789493+.1267362311985826*im;
         -.9448102919885515-.2513655659338360*im;
         -.9448102919885515+.2513655659338360*im;
         -.9044199317254038-.3718154066011398*im;
         -.9044199317254038+.3718154066011398*im;
         -.8488049242636358-.4860811780173362*im;
         -.8488049242636358+.4860811780173362*im;
         -.7788784384948876-.5922589089475263*im;
         -.7788784384948876+.5922589089475263*im;
         -.6957887129626554-.6885761217415814*im;
         -.6957887129626554+.6885761217415814*im;
         -.6009002821775153-.7734205049875106*im;
         -.6009002821775153+.7734205049875106*im;
         -.4957716929578976-.8453658760449938*im;
         -.4957716929578976+.8453658760449938*im;
         -.3821300978164688-.9031949785779086*im;
         -.3821300978164688+.9031949785779086*im;
         -.2618431760355581-.9459187037026260*im;
         -.2618431760355581+.9459187037026260*im;
         -.1368888970077769-.9727913710515751*im;
         -.1368888970077769+.9727913710515751*im;
         -.9323708922186879e-2-.9833217572983870*im;
         -.9323708922186879e-2+.9833217572983870*im;
         .1187501830426266-.9772796154416013*im;
         .1187501830426266+.9772796154416013*im;
         .2452176769219828-.9546974936920340*im;
         .2452176769219828+.9546974936920340*im;
         .3679828414167529-.9158677545985760*im;
         .3679828414167529+.9158677545985760*im;
         .4849990112751703-.8613348568446792*im;
         .4849990112751703+.8613348568446792*im;
         .5942953455083113-.7918833021200295*im;
         .5942953455083113+.7918833021200295*im;
         .6939977395882224-.7085224187149837*im;
         .6939977395882224+.7085224187149837*im;
         .7823414572835585-.6124708814099550*im;
         .7823414572835585+.6124708814099550*im;
         .8576732672059729-.5051474483453174*im;
         .8576732672059729+.5051474483453174*im;
         .9184452635324498-.3881803987843682*im;
         .9184452635324498+.3881803987843682*im;
         .9632174723875938-.2634528222232032*im;
         .9632174723875938+.2634528222232032*im;
         .9907153993859500-.1331861055266858*im;
         .9907153993859500+.1331861055266858*im;
         1.999999999999996];
    z = [z 2*ones(n)]
    p, PolyZeros(z)
end
function fl(k)
#
#  generalization of Farmer-Loizou example:
#
#    (x-1)^4k * (x-2)^3k * (x-3)^2K * (x-4)^k
#
    p = poly([ones(4*k);2*ones(3*k);3*ones(2*k);4*ones(k)])
    z = [1.0 4k; 2 4k; 3 2k; 4 k]
    p, PolyZeros(z)
end
function fl01()
#
#  generalization of Farmer-Loizou example:
#
#    (x-1)^4k * (x-2)^3k * (x-3)^2K * (X-4)^k
#
   fl(1)
end

function fl02()
#
#  generalization of Farmer-Loizou example:
#
#    (x-1)^4k * (x-2)^3k * (x-3)^2K * (X-4)^k
#
    fl(2)
end
function fl03()
#
#  generalization of Farmer-Loizou example:
#
#    (x-1)^4k * (x-2)^3k * (x-3)^2K * (X-4)^k
#
    fl(3)
end
function fl04()
#
#  generalization of Farmer-Loizou example:
#
#    (x-1)^4k * (x-2)^3k * (x-3)^2K * (X-4)^k
#
    fl(4)
end
function fl05()
#
#  generalization of Farmer-Loizou example:
#
#    (x-1)^4k * (x-2)^3k * (x-3)^2K * (X-4)^k
#
    fl(5)
end
function fl06()
#
#  generalization of Farmer-Loizou example:
#
#    (x-1)^4k * (x-2)^3k * (x-3)^2K * (X-4)^k
#
    fl(6)
end
function fl07()
#
#  generalization of Farmer-Loizou example:
#
#    (x-1)^4k * (x-2)^3k * (x-3)^2K * (X-4)^k
#
    fl(7)
end
function henrici()
#
# M. Petkovic testing polynomials, page 123
#
    p = poly([-4.1;-3.8;-2.05;-1.85;1.95;2.15;3.9;4.05])
    z = [[-4.1;-3.8;-2.05;-1.85;1.95;2.15;3.9;4.05] ones(8)];

    p, PolyZeros(z)
end
function igyp00()
#
#  generalization of Igarash and Ypma
#
    p = poly([2.35;2.37;2.39])
    z = [2.35 1;2.37 1;2.39 1];
    p, PolyZeros(z)
end
function igyp01()
#
#  generalization of Igarash and Ypma
#
    p = poly([2.35*[1;1;1];2.56])
    z = [2.35 3;2.56 1];
    p, PolyZeros(z)
end
function igyp02(m)
#
#  generalization of Igarash and Ypma
#
   p = poly([10*(1+im)*ones(m); -1*ones(10-m)])
   z = [10.0*(1+im) m; -1 10-m];
   p, PolyZeros(z)
end
function igyp02a()
#
#  generalization of Igarash and Ypma
#
    igyp02(8)
end
function igyp02b()
#
#  generalization of Igarash and Ypma
#
    igyp02(7)
end
function igyp02c()
#
#  generalization of Igarash and Ypma
#
    igyp02(3)
end
function igyp03()
#
#  generalization of Igarash and Ypma
#
    p = poly([10.0*(1+im)*[1;1;1];1;im;2;2im;3;4im;5])
    z = [10*(1+im) 3; 1 1; im 1; 2 1; 2im 1; 3 1; 4im 1; 5 1];
    p, PolyZeros(z)
end
function iliev(k)
#
#  generalization of Iliev example:
#
#    (x-1)^k (x-2)^2k (x-3)^4k
#
    p = poly([ones(k);2ones(2k);3ones(3k)])
    z = [1.0 k; 2 2k; 3 4k]
    p, PolyZeros(z)
end
function iliev01()
#
#  generalization of Iliev example:
#
#    (x-1) (x-2)^2 (x-3)^3
#
    iliev(2)
end
function iliev02()
#
#  generalization of Iliev example:
#
#    (x-1) (x-2)^2 (x-3)^3
#
    iliev(4)
end
function iliev03()
#
#  generalization of Iliev example:
#
#    (x-1) (x-2)^2 (x-3)^3
#
    iliev(8)
end
function inex(digits)
#
#  test polynomial suggested by Goedecker
#
    p = poly([(10/11)*[1;1;1;1;1]; (20/11)*[1;1;1]; (30/11)*[1;1]])
    p = round.(10^digits*reverse(p.a)) ./ 10^digits
    z = reverse([10/11 5; 20/11 3; 30/11 2]; dims=1)
    p, PolyZeros(z)
end
function inex01()
#
#  test polynomial suggested by Goedecker
#
    inex(7)
end
function inex02()
#
#  test polynomial suggested by Goedecker
#
    inex(6)
end
function inex03()
#
#  test polynomial suggested by Goedecker
#
    inex(5)
end
function inex04()
#
#  test polynomial suggested by Goedecker
#
    inex(4)
end
function jt01a()
#
#  test polynomial suggested by Jenkins and Traub
#
   a = 10.0^10
   p = [1;-1;-a^2;a^2]
   z = [a 1;-a 1;1 1]
   p, PolyZeros(z)
end
function jt01b()
#
#  test polynomial suggested by Jenkins and Traub
#
   a = 10.0^(-10)
   p = [1;-1;-a^2;a^2]
   z = [a 1;-a 1;1 1]
   p, PolyZeros(z)
end
function jt02()
#
#  test polynomial suggested by Jenkins and Traub
#
p = poly([(1.0:17)...])
    z = [(1.0:17)...]
    z = [z ones(17)]
    p, PolyZeros(z)
end
function jt03()
#
#  test polynomial suggested by Jenkins and Traub
#
    z = 1 ./ exp10.(1:8)
    p = poly(z)
    z = [z ones(8)]
    p, PolyZeros(z)
end
function jt04()
#
#  test polynomial suggested by Jenkins and Traub
#
   p = poly([0.1;0.1;0.1;0.5;0.6;0.7]);
   z = [0.5  1; 0.6  1; 0.7  1; 0.1  3];
   p, PolyZeros(z)
end
function jt05()
#
#  test polynomial suggested by Jenkins and Traub
#
    p = poly([0.1*ones(4);0.2*ones(3);0.3;0.3;0.4])
    z = [0.4  1; 0.3  2; 0.2  3; 0.1  4];
    p, PolyZeros(z)
end
function jt06()
#
#  test polynomial suggested by Jenkins and Traub
#
    p = poly([.1;1.001;.998;1.00002;.99999])
    z = [.1;1.001;.998;1.00002;.99999];
    z = [z ones(5)]
    p, PolyZeros(z)
end
function jt07a()
#
#  test polynomial suggested by Jenkins and Traub
#
   a = 10.0^(-10);
   z = [.001; .01; .1; .1+a*im; .1-a*im; 1; -10]
   p = poly(z)
   z = [z ones(7)]
   y = [-10  1; 1  1; 0.01  1; 0.001  1; 0.1  3];
   p, PolyZeros(z), PolyZeros(y)
end
function jt07b()
#
#  test polynomial suggested by Jenkins and Traub
#
   a = 1e-9
   z = [.001; .01; .1; .1+a*im; .1-a*im; 1; -10];
   p = poly(z)
   z = [z ones(7)];
   y = [-10  1; 1  1; 0.01  1; 0.001  1; 0.1  3];
   p, PolyZeros(z), PolyZeros(y)
end
function jt07c()
#
#  test polynomial suggested by Jenkins and Traub
#
    a = 1e-8
    z = [.001; .01; .1; .1+a*im; .1-a*im; 1; -10]
    p = poly(z)
    z = [z ones(7)]
    y = [-10  1; 1  1; 0.01  1; 0.001  1; 0.1  3]
    p, PolyZeros(z), PolyZeros(y)
end
function jt07d()
#
#  test polynomial suggested by Jenkins and Traub
#
    a = 1e-7
    z = [.001; .01; .1; .1+a*im; .1-a*im; 1; -10];
    p = poly(z)
    z = [z ones(7)]
    y = [-10  1; 1  1; 0.01  1; 0.001  1; 0.1  3];

    p, PolyZeros(z), PolyZeros(y)
end
function jt08()
#
#  test polynomial suggested by Jenkins and Traub
#
    p = poly([-1.0;-1;-1;-1;-1])
    z = [-1 5]
    p, PolyZeros(z)
end
function jt09()
#
#  test polynomial suggested by Jenkins and Traub
#
    p1 = [1;zeros(9);-1e-20]
    p2 = [1;zeros(9);1e20]
    p = conv(p1,p2)
    z1 = [-0.01000000000000;
          -0.00809016994375 + 0.00587785252292im;
          -0.00809016994375 - 0.00587785252292im;
          -0.00309016994375 + 0.00951056516295im;
          -0.00309016994375 - 0.00951056516295im;
          0.00309016994375 + 0.00951056516295im;
          0.00309016994375 - 0.00951056516295im;
          0.01000000000000;
          0.00809016994375 + 0.00587785252292im;
          0.00809016994375 - 0.00587785252292im;
         ]

    z2 = [-0.95105651629515 + 0.30901699437495im;
          -0.95105651629515 - 0.30901699437495im;
          -0.58778525229247 + 0.80901699437495im;
          -0.58778525229247 - 0.80901699437495im;
          -0.00000000000000 + 1.00000000000000im;
          -0.00000000000000 - 1.00000000000000im;
          0.95105651629515 + 0.30901699437495im;
          0.95105651629515 - 0.30901699437495im;
          0.58778525229247 + 0.80901699437495im;
          0.58778525229247 - 0.80901699437495im;
          ] * 100
    z = [[z1; z2] ones(20)]
    p, PolyZeros(z)
end
function jt10a()
#
#  test polynomial suggested by Jenkins and Traub
#
    a = 1e3;
    p = poly([a;1;1/a])
    z = [[a;1;1/a] ones(3)]
    p, PolyZeros(z)
end
function jt10b()
#
#  test polynomial suggested by Jenkins and Traub
#
    a = 1e6;
    p = poly([a;1;1/a])
    z = [[a;1;1/a] ones(3)]
    p, PolyZeros(z)
end
function jt10c()
#
#  test polynomial suggested by Jenkins and Traub
#
    a = 1e9;
    p = poly([a;1;1/a])
    z = [[a;1;1/a] ones(3)]
    p, PolyZeros(z)
end
function jt11a()
#
#  test polynomial suggested by Jenkins and Traub
#
    m = 15;
    y = exp.((im*pi) .* (1-m:m-1) ./ (2*m))
    z = 0.9*exp.((im*pi) .*(m:3*m) ./ (2*m))
    p1 = poly(y)
    p2 = poly(z)
    p = reverse(conv(p1.a, p2.a))
    z = [[y; z] ones(4*m)]
    p, PolyZeros(z)
end
function jt11b()
#
#  test polynomial suggested by Jenkins and Traub
#
    m = 20;
    y = exp.((im*pi) .* (1-m:m-1) ./ (2*m))
    z = 0.9*exp.((im*pi) .* (m:3*m) ./ (2*m))
    u = [y; z]
    p = poly(u)
    z = [[y; z] ones(4*m)]
    p, PolyZeros(z)
end
function jt11c()
#
#  test polynomial suggested by Jenkins and Traub
#
    m = 25;
    y = exp.((im*pi) .*  (1-m:m-1) ./ (2*m))
    z = 0.9*exp.((im*pi) .* (m:3*m) ./ (2*m))
    u = [y; z]
    p = poly(u)
    z = [u ones(length(u))]
    p, PolyZeros(z)
end
function large01()
#
#  test polynomial suggested by Jenkins and Traub
#
    p = [
         1.0;
         0;
        -0.28000000000000;
        -0.59200000000000;
        -0.96890000000000;
        -0.62020000000000;
        0.54382200000000;
        0.08842680000000;
        0.43204191000000;
        -0.15503004000000;
        0.44119337520000;
        -0.65443195750000;
        -0.35708097260000;
        -0.55597791620000;
        0.63916879370000;
        0.14637576020000;
        -0.06709193010000;
        -0.03193754010000;
        0.50857405010000;
        -0.14968258190000;
        0.63273024920000
        ];
    z = [
        -1.0 + 0.3im;
        -1.0 - 0.30im;
        -0.9 + 0.4im;
        -0.9 - 0.4im;
        -0.7 + 0.7im;
        -0.7 - 0.7im;
        -0.4 + 0.9im;
        -0.4 - 0.9im;
        -0.0 + 1.10000000000im;
        -0.000000000000 - 1.10000000000im;
        1.2;
        1.00000000000000;
        0.9 + 0.4im;
        0.9 - 0.4im;
        0.6 + 0.600000000000im;
        0.6 - 0.600000000000im;
        0.4 + 0.900000000000im;
        0.4 - 0.900000000000im;
        0.00000000000 + 0.8im
        0.00000000000 - 0.8im
        ];
    z = [z ones(20)]
    p, PolyZeros(z)
end
function large02()
#
#  test polynomial suggested by Z Zeng
#

    p, z = large01()
    p = conv(p, p)
    z.mult .*= 2
    p, z
end
function large03()
#
#  test polynomial suggested by Z Zeng
#
    p, z = large02()
    p = conv(p, p)
    z.mult .*= 2
    p, z
end
function large04()
#
#  test polynomial suggested by Z Zeng
#
    p, z = large03()
    p = conv(p, p)
    z.mult .*= 2
    p, z
end
function large05()
#
#  test polynomial suggested by Z Zeng
#
    p, z = large04()
    p = conv(p, p)
    z.mult .*= 2
    p, z
end
function largez(k::Integer)
#
#  test polynomial suggested by Z Zeng
#
    p = [
         1.0
        -0.7
        -0.19
         0.177
        -0.7364
        -0.43780
        -0.952494
        -0.2998258
        -0.00322203
        -0.328903811
        -0.4959527435
        -0.9616679762
         0.4410459281
         0.1090273141
         0.6868094008
         0.0391923826
         0.0302248540
         0.6603775863
        -0.1425784968
        -0.3437618593
         0.4357949015
        ];
    z = [
         0.5 + im
         0.5 - im
        -1 + 0.2im
        -1 - 0.2im
        -0.1 + im
        -0.1 - im
         0.8 + 0.6im
         0.8 - 0.6im
        -0.7 + 0.7im
        -0.7 - 0.7im
         1.4
        -0.4 + 0.9im
        -0.4 - 0.9im
         0.9
        -0.8 + 0.3im
        -0.8 - 0.3im
         0.3 + 0.8im
         0.3 - 0.8im
         0.6 + 0.4im
         0.6 - 0.4im
        ]

    mult = 1
    for i = 1:k
        p = conv(p, p)
        mult *= 2
    end

    z = [z mult*ones(20)]
    p, PolyZeros(z)
end

largez01() = largez(0)
largez02() = largez(1)
largez03() = largez(2)
largez04() = largez(3)
largez05() = largez(4)
largez06() = largez(5)

largez01() = largez(1)
largez02() = largez(2)
largez03() = largez(3)
largez04() = largez(4)
largez05() = largez(5)
function lgd(n)
#
#  test polynomial suggested by Goedecker
#   Legendre polynomial of degree n
#
    p0 = [1.0]
    p1 = [1.0; 0]
    if n == 0
        p = p0
    elseif n == 1
        p = p1
    elseif n > 1
        for k = 1:n-1
            q1 = (2*k+1)*[p1;0]
            m = length(q1)-length(p0)
            q0 = k*[zeros(m);p0]
            p = (q1-q0)/(k+1)
            #
            p0 = p1
            p1 = p
        end
    end

    p
end
function lgd05()
#
#  test polynomial suggested by Goedecker
#   Legendre polynomial of degree 5
#
    p = lgd(5);
    z = [0;
         0.90617984593866;
         0.53846931010568;
         -0.90617984593866;
         -0.53846931010568
        ];
    z = [z ones(5)]
    p, PolyZeros(z)
end
function lgd10()
#
#  test polynomial suggested by Goedecker
#   Legendre polynomial of degree 10
#
    p = lgd(10);
    z = [0.97390652851717
   0.86506336668898
   0.67940956829903
  -0.97390652851717
  -0.86506336668899
  -0.67940956829903
   0.43339539412925
  -0.43339539412925
   0.14887433898163
  -0.14887433898163]
    z = [z ones(10)]
    p, PolyZeros(z)
end
function lgd100()
#
#  test polynomial suggested by Goedecker
#   Legendre polynomial of degree 100
#
    p = lgd(100);
    z = [-.9997137267734412; -.9984919506395958; -.9962951347331251; -.9931249370374435;
         -.9889843952429917; -.9838775407060570; -.9778093584869183; -.9707857757637063;
         -.9628136542558155; -.9539007829254917; -.9440558701362560; -.9332885350430795;
         -.9216092981453340; -.9090295709825297; -.8955616449707270; -.8812186793850184;
         -.8660146884971646; -.8499645278795913; -.8330838798884008; -.8153892383391763;
         -.7968978923903145; -.7776279096494955; -.7575981185197072; -.7368280898020207;
         -.7153381175730564; -.6931491993558020; -.6702830156031410; -.6467619085141293;
         -.6226088602037078; -.5978474702471787; -.5725019326213812; -.5465970120650942;
         -.5201580198817631; -.4932107892081909; -.4657816497733580; -.4378974021720315;
         -.4095852916783015; -.3808729816246300; -.3517885263724217; -.3223603439005292;
         -.2926171880384720; -.2625881203715035; -.2323024818449740; -.2017898640957360;
         -.1710800805386033; -.1402031372361140; -.1091892035800611; -.7806858281343664e-1;
         -.4687168242159163e-1; -.1562898442154308e-1; .1562898442154308e-1; .4687168242159163e-1;
         .7806858281343664e-1; .1091892035800611; .1402031372361140; .1710800805386033;
         .2017898640957360; .2323024818449740; .2625881203715035; .2926171880384720;
         .3223603439005292; .3517885263724217; .3808729816246300; .4095852916783015;
         .4378974021720315; .4657816497733580; .4932107892081909; .5201580198817631;
         .5465970120650942; .5725019326213812; .5978474702471787; .6226088602037078;
         .6467619085141293; .6702830156031410; .6931491993558020; .7153381175730564;
         .7368280898020207; .7575981185197072; .7776279096494955; .7968978923903145;
         .8153892383391763; .8330838798884008; .8499645278795913; .8660146884971646;
         .8812186793850184; .8955616449707270; .9090295709825297; .9216092981453340;
         .9332885350430795; .9440558701362560; .9539007829254917; .9628136542558155;
         .9707857757637063; .9778093584869183; .9838775407060570; .9889843952429917;
         .9931249370374435; .9962951347331251; .9984919506395958; .9997137267734412]
    z = [z ones(100)]
    # Legendre polynomial, roots are sensitive, see plot.

    p, PolyZeros(z)
end
function lgd15()
#
#  test polynomial suggested by Goedecker
#   Legendre polynomial of degree 15
#
    p = lgd(15);
    z = [ 0
   0.98799251802063
   0.93727339240040
   0.84820658341064
   0.72441773136009
   0.57097217260857
   0.39415134707756
  -0.98799251802067
  -0.93727339240030
  -0.84820658341075
  -0.72441773136004
   0.20119409399743
  -0.57097217260858
  -0.39415134707756
  -0.20119409399743]
    z = [z ones(15)]

    p, PolyZeros(z)
end
function lgd20()
#
#  test polynomial suggested by Goedecker
#   Legendre polynomial of degree 20
#
    p = lgd(20);
    z = [-.9931285991850949; -.9639719272779138; -.9122344282513259; -.8391169718222188;
         -.7463319064601508; -.6360536807265150; -.5108670019508271; -.3737060887154196;
         -.2277858511416451; -.7652652113349733e-1; .7652652113349733e-1; .2277858511416451;
         .3737060887154196; .5108670019508271; .6360536807265150; .7463319064601508;
         .8391169718222188; .9122344282513259; .9639719272779138; .9931285991850949]
    z = [z ones(20)]

    p, PolyZeros(z)
end
function lgd24()
#
#  test polynomial suggested by Goedecker
#   Legendre polynomial of degree 24
#
    p = lgd(24);
    z = [-.9951872199970214; -.9747285559713095; -.9382745520027328; -.8864155270044010;
         -.8200019859739029; -.7401241915785544; -.6480936519369756; -.5454214713888395;
         -.4337935076260451; -.3150426796961634; -.1911188674736163; -.6405689286260563e-1;
         .6405689286260563e-1; .1911188674736163; .3150426796961634; .4337935076260451;
         .5454214713888395; .6480936519369756; .7401241915785544; .8200019859739029;
         .8864155270044010; .9382745520027328; .9747285559713095; .9951872199970214]
    z = [z ones(24)]
    p, PolyZeros(z)
end
function lgd50()
#
#  test polynomial suggested by Goedecker
#   Legendre polynomial of degree 50
#
    p = lgd(50);
    z = [-.9988664044200711; -.9940319694320907; -.9853540840480059; -.9728643851066921;
         -.9566109552428079; -.9366566189448779; -.9130785566557919; -.8859679795236130;
         -.8554297694299461; -.8215820708593359; -.7845558329003993; -.7444943022260685;
         -.7015524687068223; -.6558964656854394; -.6077029271849502; -.5571583045146501;
         -.5044581449074642; -.4498063349740388; -.3934143118975651; -.3355002454194374;
         -.2762881937795320; -.2160072368760418; -.1548905899981459; -.9317470156008614e-1;
         -.3109833832718888e-1; .3109833832718888e-1; .9317470156008614e-1; .1548905899981459;
         .2160072368760418; .2762881937795320; .3355002454194374; .3934143118975651;
         .4498063349740388; .5044581449074642; .5571583045146501; .6077029271849502;
         .6558964656854394; .7015524687068223; .7444943022260685; .7845558329003993;
         .8215820708593359; .8554297694299461; .8859679795236130; .9130785566557919;
         .9366566189448779; .9566109552428079; .9728643851066921; .9853540840480059;
         .9940319694320907; .9988664044200711]
    z = [z ones(50)]
    p, PolyZeros(z)
end
function miyak00()
#
#  test polynomial suggested by Goedecker
#  square of Fibocacci polynomial
#
    p = poly([(1.1+1.1im)*ones(4);(3.2+2.3im)*[1;1];2.1+1.5im]);
    z = [[(1.1+1.1im);(3.2+2.3*im);2.1+1.5*im] [4;2;1]];
    p, PolyZeros(z)
end
function miyak02()
#
#  test polynomial suggested by Goedecker
#  square of Fibocacci polynomial
#
    p = poly([(1.1+1.1im)*ones(4); (3.2+2.3im)*[1;1]; 2.1+1.5im])
    p = reverse(conv(p.a, p.a))
    z = [[(1.1+1.1im);(3.2+2.3im);2.1+1.5im] 2*[4;2;1]]
    p, PolyZeros(z)
end
function miyak04()
#
#  test polynomial suggested by Goedecker
#  square of Fibocacci polynomial
#
    p = poly([(1.1+1.1im)*ones(4); (3.2+2.3im)*[1;1]; 2.1+1.5im]);
    p = reverse(p.a)
    p = conv(p, p)
    p = conv(p, p)
    z = [[(1.1+1.1im); (3.2+2.3im); 2.1+1.5im] 4*[4;2;1]];
    p, PolyZeros(z)
end
function miyak08()
#
#  test polynomial suggested by Goedecker
#  square of Fibocacci polynomial
#
    p = poly([(1.1+1.1im)*ones(4); (3.2+2.3im)*[1;1]; 2.1+1.5im]);
    p = reverse(p.a)
    p = conv(p, p)
    p = conv(p, p)
    p = conv(p, p)
    z = [[(1.1+1.1im); (3.2+2.3im); 2.1+1.5im] 8*[4;2;1]];
    p, PolyZeros(z)
end
function near01()
#
#  test polynomial suggested by Z. Zeng
#
    e = 0.1;
    p = poly([(1-e)*ones(20); ones(20); -0.5*[1;1;1;1;1]])
    z = [1-e 20; 1 20; -0.5 5];
    p, PolyZeros(z)
end
function near02()
#
#  test polynomial suggested by Z. Zeng
#
    e = 0.01;
    p = poly([(1-e)*ones(20);ones(20);-0.5*[1;1;1;1;1]]);
    z = [1-e 20; 1 20; -0.5 5];
    p, PolyZeros(z)
end
function near03()
#
#  test polynomial suggested by Z. Zeng
#
    e = 0.001;
    p = poly([(1-e)*ones(20);ones(20);-0.5*[1;1;1;1;1]]);
    z = [1-e 20; 1 20; -0.5 5];
    p, PolyZeros(z)
end
function petk01()
#
# M. Petkovic testing polynomials (p. 109)
#
    p1 = poly([-1.0; -1; 3; 3; 3; -im*[1;1;1;1]]);
    p1 = reverse(p1.a)
    p2 = [1;-2;5];
    p2 = conv(p2,p2);
    p = conv(p1,p2);
    z = [-1 2; 3  3; -im  4; 1+2*im  2; 1-2*im  2];
    p, PolyZeros(z)
end
function petk02()
#
# M. Petkovic testing polynomials, p118
#
    p = poly([1.0;1;-im;-im;-im;5*im;5*im;-5*im;-5*im]);
    z = [1 2; -im 3; 5*im 2; -5*im 2];

    p, PolyZeros(z)
end
function petk03()
#
# M. Petkovic testing polynomials, p134
#
    p1 = [1.0;-1-99*im/70];
    p2 = [1;2;3];
    p2 = conv(p2,p2);
    p = conv(p1,p2);
    p = conv(p,[1;1]);
    z = [-1.00000000000000 + 1.41421356237309im 2;
         -1.00000000000000 - 1.41421356237309im 2;
         -1-99*im/70 1;
         -1 1];
    p, PolyZeros(z)
end
function petk04()
#
# M. Petkovic testing polynomials, p139
#
    z = [3.0; -1*[1;1;1]; 2*im*[1;1;1]; (-2+im)*[1;1]; (-2-im)*[1;1];
         (2+im)*[1;1]; (2-im)*[1;1]];
    p = poly(z);
    z = [3 1; -1 3;  2*im 3; -2+im 2; -2-im 2; 2+im 2; 2-im 2];

    p, PolyZeros(z)
end
function petk05()
#
# M. Petkovic testing polynomials, p142
#
    p = [1.0;-7;20;-28;-18;110;-92;-44;345;225];
    z = [3.00000000000000 + 0.00000000000000im
         1.00000000000000 + 2.00000000000000im
         1.00000000000000 - 2.00000000000000im
        -1.00000000000000 - 0.00000000000000im];
    z = [z [2;2;2;3]];

    p, PolyZeros(z)
end
function petk06()
#
# M. Petkovic testing polynomials, page 146
#
    y = [-1.0*[1;1;1;1];3*[1;1;1];-im;-im];

    p1 = reverse(poly(y).a)
    p2 = [1.0;-2;5];
    p2 = conv(p2,p2);
    p = conv(p1,p2);
    z = [-1.0 4; 3 3; -im 2; 1+2*im 2; 1-2*im 2];

    p, PolyZeros(z)
end
function petk07()
#
# M. Petkovic testing polynomials, page 147
#
    y = [1.0*[1;1;1];-2+im; -2-im; 5*im;5*im; -5*im; -5*im];
    p = poly(y);
    z = [1.0 3; -2+im 1;  -2-im 1;  5*im 2; -5*im 2];

    p, PolyZeros(z)
end
function toh01()
#
#  generalization of K.C. Toh and L. N. Trefethen
#
    k = [(-10:9)...]
    z = 2.0 .* (k .+ 0.5) ./ 19 + im .* sin.((2*pi) .* (k .+ 0.5) ./ 19)
    p = poly(z)
    z = [z ones(20)]

    p, PolyZeros(z)
end
function toh02()
#
#  generalization of K.C. Toh and L. N. Trefethen
#
    k = [(20:-1:0)...]
    p1 = exp10.(k)
    p2 = ones(21)
    for j = 20:-1:0, p2[21-j] = factorial(j); end
    p = p1 ./ p2;
    z = 2 .* (k .+0.5) ./ 19 + im .* sin.((2*pi) .* (k .+ 0.5) ./ 19)
    z = [1.08045842459084 + 0.92291979046775im;
        1.08045842459084 - 0.92291979046775im;
        0.57624370699834 + 0.97554880155012im;
        0.57624370699834 - 0.97554880155012im;
        0.23672740822527 + 0.94133582606706im;
        0.23672740822527 - 0.94133582606706im;
        -0.01684064074376 + 0.86388496107331im;
        -0.01684064074376 - 0.86388496107331im;
        -0.21255469312515 + 0.76040546396802im;
        -0.21255469312515 - 0.76040546396802im;
        -0.36406222857634 + 0.63986739119691im;
        -0.36406222857634 - 0.63986739119691im;
        -0.47903273599273 + 0.50773158889098im;
        -0.47903273599273 - 0.50773158889098im;
        -0.56209888019876 + 0.36770724207462im;
        -0.56209888019876 - 0.36770724207462im;
        -0.64273026408505 + 0.07449713941656im;
        -0.64273026408505 - 0.07449713941656im;
        -0.61611009709266 + 0.22255230874526im;
        -0.61611009709266 - 0.22255230874526im;
        ];
   z = [z ones(20)];
   p, PolyZeros(z)
end
function toh03()
#
#  generalization of K.C. Toh and L. N. Trefethen
#
    k = [(1:20)...]
    z = 10/11 .- exp2.(-k)
    p = poly(z)
    z = [z ones(20)]

    p, PolyZeros(z)
end
function toh04()
#
#  generalization of K.C. Toh and L. N. Trefethen
#
    p = poly((10/11)*ones(20))
    z = [10/11 20]

    p, PolyZeros(z)
end
function toh05()
#
#  generalization of K.C. Toh and L. N. Trefethen
#
    k = [(-19:0)...]
    z = exp2.(k)
    p = poly(z)
    z = [z ones(20)]

    p, PolyZeros(z)
end
function toh06a()
#
#  generalization of K.C. Toh and L. N. Trefethen
#
   p = ones(21)
   z = [0.95557280578614 + 0.29475517441090im
        0.95557280578614 - 0.29475517441090im
        0.82623877431600 + 0.56332005806362im
        0.82623877431600 - 0.56332005806362im
        0.62348980185873 + 0.78183148246803im
        0.62348980185873 - 0.78183148246803im
        0.36534102436640 + 0.93087374864420im
        0.36534102436640 - 0.93087374864420im
        0.07473009358642 + 0.99720379718118im
        0.07473009358642 - 0.99720379718118im
        -0.98883082622513 + 0.14904226617617im
        -0.98883082622513 - 0.14904226617617im
        -0.90096886790242 + 0.43388373911756im
        -0.90096886790242 - 0.43388373911756im
        -0.73305187182983 + 0.68017273777092im
        -0.73305187182983 - 0.68017273777092im
        -0.50000000000000 + 0.86602540378444im
        -0.50000000000000 - 0.86602540378444im
        -0.22252093395631 + 0.97492791218182im
        -0.22252093395631 - 0.97492791218182im];
    z = [z ones(20)]

    p, PolyZeros(z)
end
function toh06a()
#
#  generalization of K.C. Toh and L. N. Trefethen
#
    p = ones(11)
    p = conv(p,p)
    z = [0.84125353283118 + 0.54064081745560im
         0.84125353283118 - 0.54064081745560im
         0.41541501300189 + 0.90963199535452im
         0.41541501300189 - 0.90963199535452im
        -0.95949297361449 + 0.28173255684143im
        -0.95949297361449 - 0.28173255684143im
        -0.65486073394529 + 0.75574957435426im
        -0.65486073394529 - 0.75574957435426im
        -0.14231483827329 + 0.98982144188093im
        -0.14231483827329 - 0.98982144188093im];
    z = [z 2*ones(10)]

    p, PolyZeros(z)
end
function toh06c()
#
#  generalization of K.C. Toh and L. N. Trefethen
#
    p = ones(6)
    p = conv(p,p)
    p = conv(p,p);
    z = [0.50000000000000 + 0.86602540378444im
         0.50000000000000 - 0.86602540378444im
        -1.00000000000000
        -0.50000000000000 + 0.86602540378444im
        -0.50000000000000 - 0.86602540378444im]
    z = [z 4*ones(5)]

    p, PolyZeros(z)
end
function triple(k1,k2,k3)
#
#  test polynomial suggested by Goedecker
#  square of Fibocacci polynomial
#
    poly([0.9*ones(k1);ones(k2);1.1*ones(k3)])
end
function triple01()
#
#  test polynomial suggested by Goedecker
#
    p = triple(5,5,5)
    z = [[0.9;1;1.1]  [5;5;5]]

    p, PolyZeros(z)
end
function triple02()
#
#  test polynomial suggested by Goedecker
#
    p = triple(10, 10, 10)
    z = [[0.9;1;1.1]  [10; 10; 10]]

    p, PolyZeros(z)
end
function triple03()
#
#  test polynomial suggested by Goedecker
#
    p = triple(18, 10, 16)
    z = [[0.9;1;1.1]  [18; 10; 16]]

    p, PolyZeros(z)
end
function triple04()
#
#  test polynomial suggested by Goedecker
#
    p = triple(25, 15, 10)
    z = [[0.9;1;1.1]  [25; 15; 10]]

    p, PolyZeros(z)
end
function twin01()

   k = 4
   p = poly([-0.2*ones(k); 0.39*ones(k); 0.40*ones(k)])
   z = [-0.2 k; 0.39 k; 0.4 k];

   p, PolyZeros(z)
end
function twin02()
   k = 8;
   p = poly([-0.2*ones(k); 0.39*ones(k); 0.40*ones(k)])
   z = [-0.2 k; 0.39 k; 0.4 k];

   p, PolyZeros(z)
end
function twin03()

   k = 12;
   p = poly([-0.2*ones(k); 0.39*ones(k); 0.40*ones(k)])
   z = [-0.2 k; 0.39 k; 0.4 k];

   p, PolyZeros(z)
end
function twin04()

   k = 16;
   p = poly([-0.2*ones(k); 0.39*ones(k); 0.40*ones(k)])
   z = [-0.2 k; 0.39 k; 0.4 k];

   p, PolyZeros(z)
end
function uhlig01()
#
#  test polynomial used by F. Uhlig
#
    a = 0.01;
    p1 = [1;0;0;0;-a^4]
    p2 = reverse(poly([a;a;a;a]).a)
    p = conv(p1,p2)
    z = [-0.01000000000000 - 0.00000000000000im
         -0.00000000000000 + 0.01000000000000im
         -0.00000000000000 - 0.01000000000000im
         0.01000000000000 - 0.00000000000000im];
    z = [z [1;1;1;5]]
    p, PolyZeros(z)
end
function uhlig02()
#
#  test polynomial used by F. Uhlig
#
     a = 0.001;
     p1 = [1;0;0;0;-a^4]
     p2 = reverse(poly([a;a;a;a]).a)
     p = conv(p1,p2)
     z = [ -0.00000000000000 + 0.00100000000000im
         -0.00000000000000 - 0.00100000000000im
         -0.00100000000000 + 0.00000000000000im
         0.00100000000000 + 0.00000000000000im]
    z = [z [1;1;1;5]]

    p, PolyZeros(z)
end
function uhlig03()
#
#  test polynomial used by F. Uhlig
#
    r = [(3/11)*ones(12); 11/3; 11/3; (2*im/7)*ones(4); (2.5+im/4)*ones(2); 1/4]
    p = poly(r)
    z = [3/11 12; 11/3 2; 2*im/7 4; 2.5+im/4 2; 1/4 1]

    p, PolyZeros(z)
end
function uhlig04()
#
#  test polynomial used by F. Uhlig
#
    r = [(3/11)*ones(12); 11/3; 11/3; (2*im/7)*ones(4); (2.5+im/4)*ones(2); 1/8]
    p = poly(r);
    z = [3/11 12; 11/3 2; 2*im/7 4; 2.5+im/4 2; 1/8 1];

    p, PolyZeros(z)
end
function uhlig05()
#
#  test polynomial used by F. Uhlig
#
    p = poly([-1*ones(6); 2; 2]);
    z = [-1.0 6; 2 2];

    p, PolyZeros(z)
end
function uhlig06()
#
#  test polynomial used by F. Uhlig
#
    p = [1.0;0;0;0;0;0;0;0;-1];
    z = [-1.00000000000000
         -0.70710678118655 + 0.70710678118655im
         -0.70710678118655 - 0.70710678118655im
          0 + 1.00000000000000im
          0 - 1.00000000000000im
          1.00000000000000
          0.70710678118655 + 0.70710678118655im
          0.70710678118655 - 0.70710678118655im];
    z = [z ones(8)]

    p, PolyZeros(z)
end


text:
module RSA

using Random

include("utils.jl")
include("primality.jl")
include("process.jl")

struct CipherText
    data::Any
    is_str::Bool
end

"""
RSA Crypto RSASystem
"""
struct RSASystem{T <: Integer}
    bit_size::T
    primality_test::Function
    p::BigInt
    q::BigInt
    n::BigInt
    phi::BigInt
    public_key::Tuple{BigInt, BigInt}
    private_key::Tuple{BigInt, BigInt}
end

function RSASystem(bit_size::T = 256, primality_test::Function = miller_rabin) where {T <: Integer}
    p = generate_large_prime(bit_size, primality_test)
    q = generate_large_prime(bit_size, primality_test)

    while p == q
        q = generate_large_prime(bit_size, primality_test)
    end
    n = p * q
    phi = (p - 1) * (q - 1)
    public_key, private_key = generate_keys(n, phi)
    RSASystem(bit_size, primality_test, p, q, n, phi, public_key, private_key)
end

"""
Generates public and private keys
"""
function generate_keys(n::BigInt, phi::BigInt)
    e = rand(2:phi-1)

    while gcd(e, phi) != 1
        e = rand(2:phi-1)
    end

    d = inverse(e, phi)
    public_key = (e, n)
    private_key = (d, n)
    public_key, private_key
end

include("encrypt.jl")

export encrypt!, decrypt!, RSASystem

end # module


text:
# Word Problems Logical Reasoning Tests
# http://www.graduatewings.co.uk/word-problems-logical-tests
using AbstractLogic

# Question 1
# Sam is shorter than Hannah. William is shorter than Michael who is shorter than
# Hannah. Hannah and William are shorter than Fiona. Who is the tallest?
logicset = logicalparse(["Sam, Hannah, William, Michael, Fiona ∈ 1:5"])
logicset = logicalparse(["Sam < Hannah; William < Michael; Michael < Hannah"], logicset)
logicset = logicalparse(["Hannah < Fiona; William < Fiona"], logicset)

# A) Sam    B) Hannah    C) William    D) Fiona    E) Cannot Say
search("{{i}} > {{!i}}", logicset)
# :Fiona is a match with 13 feasible combinations out of 13.
# Therefore we choose D) Fiona

# Question 2
# Ben is taller than Sharon who is taller than Ralph. Ben is also taller than Mike.
logicset = logicalparse(["Ben, Sharon, Ralph, Mike ∈ 1:4"])
logicset = logicalparse(["Ben > Sharon; Sharon > Ralph; Ben > Mike"], logicset)

# Who is the shortest?
# A) Ben    B) Sharon    C) Ralph    D) Mike    E) Cannot Say
search("{{i}} < {{!i}}", logicset)
# :Ralph is a possible match with 6 feasible combinations out of 11.
# :Mike is a possible match with 1 feasible combinations out of 11.
# Therefore we choose E) Cannot Say


text:
module StarGene

# package code goes here

include("stargene.jl")
include("samplegene.jl")

end # module


text:
abstract type AbstractSkyDomain end

struct LowerHemisphere <: AbstractSkyDomain end
struct BothHemispheres <: AbstractSkyDomain end

abstract type AbstractDirectionSampler{AbstractSkyDomain} end

struct RandomSampler{D} <: AbstractDirectionSampler{D}
    RandomSampler(domain::AbstractSkyDomain) = new{domain}()
    RandomSampler(; domain::AbstractSkyDomain = LowerHemisphere()) = new{typeof(domain)}()
end
struct EvenSampler{D} <: AbstractDirectionSampler{D}
    EvenSampler(domain::AbstractSkyDomain) = new{domain}()
    EvenSampler(; domain::AbstractSkyDomain = LowerHemisphere()) = new{typeof(domain)}()
end
struct WeierstrassSampler{D} <: AbstractDirectionSampler{D}
    resolution::Float64
    WeierstrassSampler(res, domain::AbstractSkyDomain) = new{typeof(domain)}(res)
    WeierstrassSampler(; res = 100.0, domain::AbstractSkyDomain = LowerHemisphere()) =
        WeierstrassSampler(res, domain)
end

sample_azimuth(sm::AbstractDirectionSampler{D}, i) where {D} =
    error("Not implemented for $(typeof(sm)).")
sample_elevation(sm::AbstractDirectionSampler{D}, i) where {D} =
    error("Not implemented for $(typeof(sm)).")
sample_angles(sm::AbstractDirectionSampler{D}, i, N) where {D} =
    (sample_elevation(sm, i / N), sample_azimuth(sm, i))

# even
sample_azimuth(::RandomSampler{D}, i) where {D} = rand(Float64) * 2π
sample_elevation(::RandomSampler{LowerHemisphere}, i) = acos(1 - rand(Float64))
sample_elevation(::RandomSampler{BothHemispheres}, i) = acos(1 - 2 * rand(Float64))

# even
sample_azimuth(::EvenSampler{D}, i) where {D} = π * (1 + √5) * i
sample_elevation(::EvenSampler{LowerHemisphere}, i) = acos(1 - i)
sample_elevation(::EvenSampler{BothHemispheres}, i) = acos(1 - 2i)

# uniform
sample_azimuth(::WeierstrassSampler{D}, i) where {D} = π * (1 + √5) * i
sample_elevation(sm::WeierstrassSampler{LowerHemisphere}, i) =
    π - 2atan(√(sm.resolution / i))
function sample_elevation(sm::WeierstrassSampler{BothHemispheres}, i)
    ϕ = 2atan(√(sm.resolution / i))
    if iseven(i)
        ϕ
    else
        π - ϕ
    end
end
sample_angles(sm::WeierstrassSampler{D}, i, N) where {D} =
    (sample_elevation(sm, i), sample_azimuth(sm, i))

function convert_angles(a, r, θ, ϕ, θ_obs, ϕ_obs)
    ϕ̃ = ϕ - ϕ_obs
    R = sqrt(r^2 + a^2)
    o1 = r * R * sin(θ) * sin(θ_obs) * cos(ϕ̃) + R^2 * cos(θ) * cos(θ_obs)
    o2 = R * cos(θ) * sin(θ_obs) * cos(ϕ̃) - r * sin(θ) * cos(θ_obs)
    o3 = sin(θ_obs) * sin(ϕ̃) / sin(θ)
    sigma = r^2 + a^2 * cos(θ)^2
    -o1 / sigma, -o2 / sigma, o3 / R
end

function vector_to_local_sky(m::AbstractMetricParams{T}, u, θ, ϕ) where {T}
    error("Not implemented for $(typeof(m))")
end

function vector_to_local_sky(m::BoyerLindquist{T}, u, θ, ϕ) where {T}
    convert_angles(m.a, u[2], u[3], u[4], θ, ϕ)
end

export LowerHemisphere, BothHemispheres, RandomSampler, EvenSampler, WeierstrassSampler


text:
using Flux, Zygote
using Revise
using MFGnet
using Test

R = Float32
m = 16
d = 2
nTh = 2

(K1,b1) = (R(.1)*randn(R,m,d+1),R(.1)*randn(R,m))
Θ1 = (K1,b1)
(K2,b2) = (R(.1)*repeat(randn(R,m,m),1,1,nTh),R(.1)*repeat(randn(R,m),1,nTh))
Θ2 = (K2,b2)
ΘN = (Θ1, Θ2)
(w0,A0,c0) = (ones(R,m)/R(sqrt(m)),zeros(R,d+1,d+1),zeros(R,d+1))
Θ = (w0,ΘN,A0,c0)

As = param(copy(A0))
K1s = param(copy(K1))
K2s = param(copy(K2))
b1s = param(copy(b1))
b2s = param(copy(b2))
ws = param(copy(w0))
cs = param(copy(c0))

parms = (ws,((K1s,b1s),(K2s,b2s)),As,cs)
Θopt  = Array(MFGnet.param2vec(parms))

@test MFGnet.getParmsType(parms)==R
@test typeof(Θopt[1])==R
# println("Param2Vec and getParmsType Test Passed!")
#


text:
@testset "trainable" begin
    ν = 2.0; c = 3.0; d = 2.0; γ = 2.0; α = 2.5; h = 0.5; r = rand(3)

    function test_params(kernel, reference)
        params_kernel = params(kernel)
        params_reference = params(reference)

        @test length(params_kernel) == length(params_reference)
        @test all(p == q for (p, q) in zip(params_kernel, params_reference))
    end

    kc = ConstantKernel(c=c)
    test_params(kc, ([c],))

    kfbm = FBMKernel(h = h)
    test_params(kfbm, ([h],))

    kge = GammaExponentialKernel(γ=γ)
    test_params(kge, ([γ],))

    kgr = GammaRationalQuadraticKernel(γ=γ, α=α)
    test_params(kgr, ([α], [γ]))

    kl = LinearKernel(c=c)
    test_params(kl, ([c],))

    km = MaternKernel(ν=ν)
    test_params(km, ([ν],))

    kp = PolynomialKernel(c=c, d=d)
    test_params(kp, ([d], [c]))

    kpe = PeriodicKernel(r = r)
    test_params(kpe, (r,))

    kr = RationalQuadraticKernel(α=α)
    test_params(kr, ([α],))

    k = km + kc
    test_params(k, (k.weights, km, kc))

    k = km * kc
    test_params(k, (km, kc))

    s = 2.0
    k = transform(km, s)
    test_params(k, ([s], km))

    v = [2.0]
    k = transform(kc, v)
    test_params(k, (v, kc))

    P = rand(3, 2)
    k = transform(km, LinearTransform(P))
    test_params(k, (P, km))

    k = transform(km, LinearTransform(P) ∘ ScaleTransform(s))
    test_params(k, ([s], P, km))

    c = Chain(Dense(3, 2))
    k = transform(km, FunctionTransform(c))
    test_params(k, (c, km))
end


text:
# Display an array in the form
# [ -3.0e-01 -5.1e-01  1.9e-01 ... -2.3e-01 -4.4e-01  2.4e-01 ]
# with ndisp/2 elements on each side.
function vec2str(x :: Array{Float64,1}; ndisp :: Int=7)
  n = length(x);
  if n <= ndisp
    ndisp = n;
    nside = n;
  else
    nside = max(1, div(ndisp - 1, 2));
  end
  s = "[ ";
  i = 1;
  while i <= nside
    s *= @sprintf("%8.1e ", x[i]);
    i += 1;
  end
  if i <= div(n, 2)
    s *= "... ";
  end;
  i = max(i, n - nside + 1);
  while i <= n
    s *= @sprintf("%8.1e ", x[i]);
    i += 1;
  end
  s *= "]";
  return s;
end


text:
module MakeLatexLambda

using JSON
using AWSCore
using AWSS3
# using AWSLambda: create_lambda, update_lambda, lambda_configuration,
#                 invoke_lambda
using InfoZIP

function all()
    download()
    build()
    zip()
    deploy()
end

# Download TeXLive installer.
function download()
    run(`wget
         http://mirror.ctan.org/systems/texlive/tlnet/install-tl-unx.tar.gz`)
    print("Download completed")
end

# Build docker image from Dockerfile.
function build()
    run(`sudo docker build -t octech/lambdalatex .`)
    print("Build completed")
end

# Create Lambda deployment .ZIP file from docker image.
function zip()
    rm("latexlambda.zip", force=true)
    run(`docker run --rm -it -v $(pwd()):/var/host octech/lambdalatex zip --symlinks -r -9 /var/host/latexlambda.zip .`)
end

# Deploy .ZIP file to S3 bucket
function deploy()
    AWSCore.set_debug_level(2)
    aws = AWSCore.default_aws_config()
    n = AWSCore.aws_account_number(aws)
    aws[:lambda_bucket] = "tuftelatex"
    s3_create_bucket(aws[:lambda_bucket])
    s3_put(aws[:lambda_bucket], "latexlambda.zip", read("latexlambda.zip"))
end

# Docker image interactive shell.
function shell()
    run(`docker run --rm -it -v $(pwd()):/var/host octech/lambdalatex bash`)
end

test_zip = base64encode(create_zip("document.tex" =>
                                   readstring("test_input.tex")))

# Test latex in local docker image.
function localtest()
    write("test_input.zip.base64", test_zip)
    pycmd = """
        import lambda_function
        import json
        import base64
        out = lambda_function.lambda_handler(
              {'input': open('/var/host/test_input.zip.base64').read()}, {})
        with open('/var/host/test_output.json', 'w') as f:
            f.write(json.dumps(out))
        """
    run(`docker run --security-opt seccomp:unconfined --rm -v $(pwd()):/var/host octech/lambdalatex strace -f -t -e trace=file -o /var/host/strace.out python3 -c $pycmd`)
    out = JSON.parse(readstring("test_output.json"))
    write("test_output_local.pdf", base64decode(out["output"]))
    write("test_output_local.stdout", out["stdout"])
end

# Remove intermediate files.
function clean()
    rm("test_output.json", force=true)
    rm("test_output_local.stdout", force=true)
    rm("test_output_local.pdf", force=true)
    rm("test_output_lambda.stdout", force=true)
    rm("test_output_lambda.pdf", force=true)
    rm("latexlambda.zip", force=true)
end

if length(ARGS) == 0
    all()
else
    include_string("$(ARGS[1])()")
end

end # module MakeLatexLambda

max_stars_repo_name: PrimerLabs/tuftelambdalatex

max_issues_repo_name: PrimerLabs/tuftelambdalatex

max_forks_repo_name: PrimerLabs/tuftelambdalatex