text:
# Compile-time constants
# From HPE MPT 2.23's `mpi.h` and `mpio.h` on HLRS Hawk (2022-04-22)

# Implementation limits:
const MPI_MAX_DATAREP_STRING = Cint(128)
const MPI_MAX_ERROR_STRING = Cint(256)
const MPI_MAX_INFO_KEY = Cint(255)
const MPI_MAX_INFO_VAL = Cint(1024)
const MPI_MAX_LIBRARY_VERSION_STRING = Cint(256)
const MPI_MAX_OBJECT_NAME = Cint(128)
const MPI_MAX_PORT_NAME = Cint(256)
const MPI_MAX_PROCESSOR_NAME = Cint(256)

# Types

# Various (signed) integer types:
const MPI_Aint   = Clong
const MPI_Fint   = Cint
const MPI_Count  = Clonglong
const MPI_Offset = Clonglong

# Status:
struct MPI_Status
    MPI_SOURCE::Cint
    MPI_TAG::Cint
    MPI_ERROR::Cint
    _private0::Cint
    _private1::Clong
end

# Opaque handles:
const MPI_Comm       = Cuint
const MPI_Datatype   = Cuint
const MPI_Errhandler = Cuint
const MPI_File       = Ptr{Cvoid}
const MPI_Group      = Cuint
const MPI_Info       = Cuint
const MPI_Message    = Cint
const MPI_Op         = Cuint
const MPI_Request    = Cuint
const MPI_Win        = Cuint

# Function pointers:
const MPI_Comm_copy_attr_function = Ptr{Cvoid}
const MPI_Comm_delete_attr_function = Ptr{Cvoid}
const MPI_Comm_errhandler_function = Ptr{Cvoid}
const MPI_Comm_errhandler_fn = MPI_Comm_errhandler_function
const MPI_Copy_function = Ptr{Cvoid}
const MPI_Datarep_conversion_function = Ptr{Cvoid}
const MPI_Datarep_extent_function = Ptr{Cvoid}
const MPI_Delete_function = Ptr{Cvoid}
const MPI_File_errhandler_function = Ptr{Cvoid}
const MPI_File_errhandler_fn = MPI_File_errhandler_function
const MPI_Grequest_cancel_function = Ptr{Cvoid}
const MPI_Grequest_free_function = Ptr{Cvoid}
const MPI_Grequest_query_function = Ptr{Cvoid}
const MPI_Type_copy_attr_function = Ptr{Cvoid}
const MPI_Type_delete_attr_function = Ptr{Cvoid}
const MPI_User_function = Ptr{Cvoid}
const MPI_Win_copy_attr_function = Ptr{Cvoid}
const MPI_Win_delete_attr_function = Ptr{Cvoid}
const MPI_Win_errhandler_function = Ptr{Cvoid}
const MPI_Win_errhandler_fn = MPI_Win_errhandler_function

# Run-time constants

@const_ref MPI_ANY_SOURCE Cint -2
@const_ref MPI_ANY_TAG    Cint -1
@const_ref MPI_PROC_NULL  Cint -1
@const_ref MPI_ROOT       Cint -4

@const_ref MPI_GRAPH      Cint 1
@const_ref MPI_CART       Cint 2
@const_ref MPI_DIST_GRAPH Cint 3

@const_ref MPI_IDENT      Cint 0
@const_ref MPI_CONGRUENT  Cint 1
@const_ref MPI_SIMILAR    Cint 2
@const_ref MPI_UNEQUAL    Cint 3

@const_ref MPI_KEYVAL_INVALID     Cint  0
@const_ref MPI_UNDEFINED          Cint -3

@const_ref MPI_TAG_UB             Cint  1
@const_ref MPI_HOST               Cint  2
@const_ref MPI_IO                 Cint  3
@const_ref MPI_WTIME_IS_GLOBAL    Cint  4
@const_ref MPI_UNIVERSE_SIZE      Cint  9
@const_ref MPI_LASTUSEDCODE       Cint 13
@const_ref MPI_APPNUM             Cint 11

@const_ref MPI_WIN_BASE           Cint 1
@const_ref MPI_WIN_SIZE           Cint 2
@const_ref MPI_WIN_DISP_UNIT      Cint 3
@const_ref MPI_WIN_CREATE_FLAVOR  Cint 7
@const_ref MPI_WIN_MODEL          Cint 9

@const_ref MPI_COMBINER_NAMED            Cint -1
@const_ref MPI_COMBINER_DUP              Cint  7
@const_ref MPI_COMBINER_CONTIGUOUS       Cint  0
@const_ref MPI_COMBINER_VECTOR           Cint  1
@const_ref MPI_COMBINER_HVECTOR_INTEGER  Cint 12
@const_ref MPI_COMBINER_HVECTOR          Cint  2
@const_ref MPI_COMBINER_INDEXED          Cint  3
@const_ref MPI_COMBINER_HINDEXED_INTEGER Cint 11
@const_ref MPI_COMBINER_HINDEXED         Cint  4
@const_ref MPI_COMBINER_INDEXED_BLOCK    Cint 13
@const_ref MPI_COMBINER_STRUCT_INTEGER   Cint 15
@const_ref MPI_COMBINER_STRUCT           Cint  5
@const_ref MPI_COMBINER_SUBARRAY         Cint 16
@const_ref MPI_COMBINER_DARRAY           Cint  6
@const_ref MPI_COMBINER_F90_REAL         Cint 10
@const_ref MPI_COMBINER_F90_COMPLEX      Cint  8
@const_ref MPI_COMBINER_F90_INTEGER      Cint  9
@const_ref MPI_COMBINER_RESIZED          Cint 14
@const_ref MPI_COMBINER_HINDEXED_BLOCK   Cint 17

@const_ref MPI_COMM_TYPE_SHARED Cint 1

@const_ref MPI_DISTRIBUTE_BLOCK  Cint 121
@const_ref MPI_DISTRIBUTE_CYCLIC Cint 122
@const_ref MPI_DISTRIBUTE_NONE   Cint 123

@const_ref MPI_SUCCESS                   Cint  0
@const_ref MPI_ERR_BUFFER                Cint  1
@const_ref MPI_ERR_COUNT                 Cint  2
@const_ref MPI_ERR_TYPE                  Cint  3
@const_ref MPI_ERR_TAG                   Cint  4
@const_ref MPI_ERR_COMM                  Cint  5
@const_ref MPI_ERR_RANK                  Cint  6
@const_ref MPI_ERR_ROOT                  Cint  8
@const_ref MPI_ERR_TRUNCATE              Cint 15
@const_ref MPI_ERR_GROUP                 Cint  9
@const_ref MPI_ERR_OP                    Cint 10
@const_ref MPI_ERR_REQUEST               Cint  7
@const_ref MPI_ERR_TOPOLOGY              Cint 11
@const_ref MPI_ERR_DIMS                  Cint 12
@const_ref MPI_ERR_ARG                   Cint 13
@const_ref MPI_ERR_OTHER                 Cint 16
@const_ref MPI_ERR_UNKNOWN               Cint 14
@const_ref MPI_ERR_INTERN                Cint 17
@const_ref MPI_ERR_IN_STATUS             Cint 18
@const_ref MPI_ERR_PENDING               Cint 19
@const_ref MPI_ERR_ACCESS                Cint 28
@const_ref MPI_ERR_AMODE                 Cint 29
@const_ref MPI_ERR_BAD_FILE              Cint 31
@const_ref MPI_ERR_CONVERSION            Cint 33
@const_ref MPI_ERR_DUP_DATAREP           Cint 35
@const_ref MPI_ERR_FILE_EXISTS           Cint 36
@const_ref MPI_ERR_FILE_IN_USE           Cint 37
@const_ref MPI_ERR_FILE                  Cint 38
@const_ref MPI_ERR_IO                    Cint 43
@const_ref MPI_ERR_NO_SPACE              Cint 49
@const_ref MPI_ERR_NO_SUCH_FILE          Cint 50
@const_ref MPI_ERR_READ_ONLY             Cint 53
@const_ref MPI_ERR_UNSUPPORTED_DATAREP   Cint 59
@const_ref MPI_ERR_INFO                  Cint 42
@const_ref MPI_ERR_INFO_KEY              Cint 39
@const_ref MPI_ERR_INFO_VALUE            Cint 41
@const_ref MPI_ERR_INFO_NOKEY            Cint 40
@const_ref MPI_ERR_NAME                  Cint 46
@const_ref MPI_ERR_NO_MEM                Cint 47
@const_ref MPI_ERR_NOT_SAME              Cint 48
@const_ref MPI_ERR_PORT                  Cint 51
@const_ref MPI_ERR_QUOTA                 Cint 52
@const_ref MPI_ERR_SERVICE               Cint 56
@const_ref MPI_ERR_SPAWN                 Cint 58
@const_ref MPI_ERR_UNSUPPORTED_OPERATION Cint 60
@const_ref MPI_ERR_WIN                   Cint 61
@const_ref MPI_ERR_BASE                  Cint 32
@const_ref MPI_ERR_LOCKTYPE              Cint 45
@const_ref MPI_ERR_KEYVAL                Cint 44
@const_ref MPI_ERR_RMA_CONFLICT          Cint 54
@const_ref MPI_ERR_RMA_SYNC              Cint 55
@const_ref MPI_ERR_SIZE                  Cint 57
@const_ref MPI_ERR_DISP                  Cint 34
@const_ref MPI_ERR_ASSERT                Cint 30
@const_ref MPI_ERR_RMA_RANGE             Cint 62
@const_ref MPI_ERR_RMA_ATTACH            Cint 63
@const_ref MPI_ERR_RMA_SHARED            Cint 64
@const_ref MPI_ERR_RMA_FLAVOR            Cint 65

@const_ref MPI_LOCK_EXCLUSIVE Cint 2
@const_ref MPI_LOCK_SHARED    Cint 1

@const_ref MPI_MODE_RDONLY           Cint   2
@const_ref MPI_MODE_RDWR             Cint   8
@const_ref MPI_MODE_WRONLY           Cint   4
@const_ref MPI_MODE_CREATE           Cint   1
@const_ref MPI_MODE_EXCL             Cint  64
@const_ref MPI_MODE_DELETE_ON_CLOSE  Cint  16
@const_ref MPI_MODE_UNIQUE_OPEN      Cint  32
@const_ref MPI_MODE_APPEND           Cint 128
@const_ref MPI_MODE_SEQUENTIAL       Cint 256

@const_ref MPI_MODE_NOCHECK          Cint  1
@const_ref MPI_MODE_NOSTORE          Cint  2
@const_ref MPI_MODE_NOPUT            Cint  4
@const_ref MPI_MODE_NOPRECEDE        Cint  8
@const_ref MPI_MODE_NOSUCCEED        Cint 16

@const_ref MPI_ORDER_C       Cint 56
@const_ref MPI_ORDER_FORTRAN Cint 57

@const_ref MPI_SEEK_SET Cint 600
@const_ref MPI_SEEK_CUR Cint 602
@const_ref MPI_SEEK_END Cint 604

@const_ref MPI_THREAD_SINGLE     Cint 0
@const_ref MPI_THREAD_FUNNELED   Cint 1
@const_ref MPI_THREAD_SERIALIZED Cint 2
@const_ref MPI_THREAD_MULTIPLE   Cint 3

@const_ref MPI_TYPECLASS_REAL    Cint 2
@const_ref MPI_TYPECLASS_INTEGER Cint 1
@const_ref MPI_TYPECLASS_COMPLEX Cint 3

@const_ref MPI_ARGV_NULL      Ptr{Cvoid} C_NULL
@const_ref MPI_ARGVS_NULL     Ptr{Cvoid} C_NULL

@const_ref MPI_UNWEIGHTED     Ptr{Cvoid} cglobal((:MPI_UNWEIGHTED, libmpi), Ptr{Cvoid})
@const_ref MPI_WEIGHTS_EMPTY  Ptr{Cvoid} cglobal((:MPI_WEIGHTS_EMPTY, libmpi), Ptr{Cvoid})
@const_ref MPI_BOTTOM         Ptr{Cvoid} C_NULL
@const_ref MPI_IN_PLACE       Ptr{Cvoid} cglobal((:mpi_sgi_inplace, libmpi))

@const_ref MPI_COMM_NULL  MPI_Comm 0
@const_ref MPI_COMM_SELF  MPI_Comm 2
@const_ref MPI_COMM_WORLD MPI_Comm 1

@const_ref MPI_COMM_DUP_FN           MPI_Comm_copy_attr_function   cglobal((:MPI_COMM_DUP_FN, libmpi))
@const_ref MPI_COMM_NULL_COPY_FN     MPI_Comm_copy_attr_function   cglobal((:MPI_COMM_NULL_COPY_FN, libmpi))
@const_ref MPI_COMM_NULL_DELETE_FN   MPI_Comm_delete_attr_function cglobal((:MPI_COMM_NULL_DELETE_FN, libmpi))

@const_ref MPI_DATATYPE_NULL MPI_Datatype 0

# Only define C constants, as we don't need the Fortran or C++ ones (and Fortran ones are compiler-dependent).
@const_ref MPI_CHAR                MPI_Datatype  1
@const_ref MPI_UNSIGNED_CHAR       MPI_Datatype  5
@const_ref MPI_SHORT               MPI_Datatype  2
@const_ref MPI_UNSIGNED_SHORT      MPI_Datatype  6
@const_ref MPI_INT                 MPI_Datatype  3
@const_ref MPI_UNSIGNED            MPI_Datatype  7
@const_ref MPI_LONG                MPI_Datatype  4
@const_ref MPI_UNSIGNED_LONG       MPI_Datatype  8
# @const_ref MPI_LONG_LONG_INT       MPI_Datatype
@const_ref MPI_FLOAT               MPI_Datatype  9
@const_ref MPI_DOUBLE              MPI_Datatype 10
@const_ref MPI_LONG_DOUBLE         MPI_Datatype 11
@const_ref MPI_LONG_LONG           MPI_Datatype 12
@const_ref MPI_BYTE                MPI_Datatype 27
@const_ref MPI_WCHAR               MPI_Datatype 40
@const_ref MPI_PACKED              MPI_Datatype 28

@const_ref MPI_LB                  MPI_Datatype 30
@const_ref MPI_UB                  MPI_Datatype 29

@const_ref MPI_FLOAT_INT           MPI_Datatype 31
@const_ref MPI_DOUBLE_INT          MPI_Datatype 32
@const_ref MPI_LONG_INT            MPI_Datatype 33
@const_ref MPI_SHORT_INT           MPI_Datatype 35
@const_ref MPI_LONG_DOUBLE_INT     MPI_Datatype 36

@const_ref MPI_2INT                MPI_Datatype 34

@const_ref MPI_SIGNED_CHAR         MPI_Datatype 41
@const_ref MPI_UNSIGNED_LONG_LONG  MPI_Datatype 42

@const_ref MPI_INT8_T                   MPI_Datatype 47
@const_ref MPI_INT16_T                  MPI_Datatype 48
@const_ref MPI_INT32_T                  MPI_Datatype 49
@const_ref MPI_INT64_T                  MPI_Datatype 50
@const_ref MPI_UINT8_T                  MPI_Datatype 51
@const_ref MPI_UINT16_T                 MPI_Datatype 52
@const_ref MPI_UINT32_T                 MPI_Datatype 53
@const_ref MPI_UINT64_T                 MPI_Datatype 54
@const_ref MPI_C_BOOL                   MPI_Datatype 57
@const_ref MPI_C_FLOAT_COMPLEX          MPI_Datatype 58
@const_ref MPI_C_DOUBLE_COMPLEX         MPI_Datatype 59
@const_ref MPI_C_LONG_DOUBLE_COMPLEX    MPI_Datatype 60
@const_ref MPI_AINT                     MPI_Datatype 55
@const_ref MPI_OFFSET                   MPI_Datatype 56
@const_ref MPI_COUNT                    MPI_Datatype 61
# @const_ref MPIX_C_FLOAT16               MPI_Datatype

# aliases
const MPI_LONG_LONG_INT = MPI_LONG_LONG
const MPI_C_COMPLEX = MPI_C_FLOAT_COMPLEX

@const_ref MPI_ERRHANDLER_NULL  MPI_Errhandler 0
@const_ref MPI_ERRORS_ARE_FATAL MPI_Errhandler 1
@const_ref MPI_ERRORS_RETURN    MPI_Errhandler 2
# @const_ref MPI_ERRORS_ABORT MPI_Errhandler

@const_ref MPI_FILE_NULL MPI_File C_NULL

@const_ref MPI_GROUP_EMPTY MPI_Group 1
@const_ref MPI_GROUP_NULL  MPI_Group 0

@const_ref MPI_INFO_ENV  MPI_Info 1
@const_ref MPI_INFO_NULL MPI_Info 0

@const_ref MPI_MESSAGE_NO_PROC MPI_Message -1
@const_ref MPI_MESSAGE_NULL    MPI_Message  0

@const_ref MPI_DISPLACEMENT_CURRENT MPI_Offset -54278278

@const_ref MPI_OP_NULL MPI_Op  0
@const_ref MPI_MAX     MPI_Op  1
@const_ref MPI_MIN     MPI_Op  2
@const_ref MPI_SUM     MPI_Op  3
@const_ref MPI_PROD    MPI_Op  4
@const_ref MPI_LAND    MPI_Op  5
@const_ref MPI_BAND    MPI_Op  6
@const_ref MPI_LOR     MPI_Op  7
@const_ref MPI_BOR     MPI_Op  8
@const_ref MPI_LXOR    MPI_Op  9
@const_ref MPI_BXOR    MPI_Op 10
@const_ref MPI_MINLOC  MPI_Op 12
@const_ref MPI_MAXLOC  MPI_Op 11
@const_ref MPI_REPLACE MPI_Op 13
@const_ref MPI_NO_OP   MPI_Op 14

@const_ref MPI_REQUEST_NULL MPI_Request 0

@const_ref MPI_STATUS_IGNORE    Ptr{Cvoid} cglobal((:mpi_sgi_status_ignore, libmpi))
@const_ref MPI_STATUSES_IGNORE  Ptr{Cvoid} cglobal((:mpi_sgi_status_ignore, libmpi))

@const_ref MPI_TYPE_DUP_FN          MPI_Type_copy_attr_function   cglobal((:MPI_TYPE_DUP_FN, libmpi))
@const_ref MPI_TYPE_NULL_COPY_FN    MPI_Type_copy_attr_function   cglobal((:MPI_TYPE_NULL_COPY_FN, libmpi))
@const_ref MPI_TYPE_NULL_DELETE_FN  MPI_Type_delete_attr_function cglobal((:MPI_TYPE_NULL_DELETE_FN, libmpi))

@const_ref MPI_WIN_NULL MPI_Win 0

@const_ref MPI_WIN_DUP_FN           MPI_Win_copy_attr_function   cglobal((:MPI_WIN_DUP_FN, libmpi))
@const_ref MPI_WIN_NULL_COPY_FN     MPI_Win_copy_attr_function   cglobal((:MPI_WIN_NULL_COPY_FN, libmpi))
@const_ref MPI_WIN_NULL_DELETE_FN   MPI_Win_delete_attr_function cglobal((:MPI_WIN_NULL_DELETE_FN, libmpi))


text:
#!/usr/bin/env julia

__precompile__(true)
module Fxcrd

include("./ShellInterop.jl")
include("./Tokenize.jl")
include("./ShellCompiler.jl")
include("./MacroDatabase.jl")
include("./MacroExpander.jl")
include("./StartMacros.jl")
include("./Source.jl")
include("./Interactive.jl")

using .Tokenize
using .MacroExpander
using .ShellInterop
using .Source

function main(ev,rd,files...)
  try
    configure_shell(ev,rd)
    source(files...)
  catch err
    shell(err_sh)
    shell(end_sh)
    rethrow(err)
    exit()
  finally
    shell(end_sh)
  end
end

greenprint(x) = print_with_color(:green, x, "\n")

function interactive()
  println()
  greenprint("  ------------------")
  greenprint("   --- hello !! --- ")
  greenprint("  ------------------")
  println()

  open(`fxcrd_compile --print --send`) do compiler
    ev = readline(compiler, chomp=true)
    rd = readline(compiler, chomp=true)
    configure_shell(ev, rd)

    shell() do sh
      eval_sh(sh, """
        hi()    { echo "\033[32;1m :)\033[0m"; }
        hello() { echo "\033[32;1m :)\033[0m"; }
      """)
      Interactive.main_loop(sh)
      greenprint(" bye bye! ")
      println()
      end_sh(sh)
    end
  end
end

end

ext: jl

max_stars_repo_name: fcard/Fxcrd

max_issues_repo_name: fcard/Fxcrd

max_forks_repo_name: fcard/Fxcrd


text:
using DelimitedFiles
using Random
using Statistics

function sigmoid(x)
        return 1/(1+exp(-x));
end

function sigmoidDerivative(x)
        return sigmoid(x)*(1 - sigmoid(x));
end

data = readdlm(ARGS[1], ',', Float64, '\n');
numOfDatapoints = size(data,2)
hiddenLayerSize = 15;
learningRate = 0.2;

resistanceMean = mean(data[1,:]);
resistanceStdDev = std(data[1,:]);
capacitanceMean = mean(data[2,:]);
capacitanceStdDev = std(data[2,:]);
timeMean = mean(data[3,:]);
timeStdDev = std(data[3,:]);

inputs = vcat(data[1:3,:], ones(numOfDatapoints)'); # Adds row of ones for bias term in weights
expectedOutput = data[4,:];

# Normalize inputs
inputs[1,:] .-= resistanceMean;
inputs[1,:] ./= resistanceStdDev;
inputs[2,:] .-= capacitanceMean;
inputs[2,:] ./= capacitanceStdDev;
inputs[3,:] .-= timeMean;
inputs[3,:] ./= timeStdDev;

global layerOneWeights = rand!(zeros(hiddenLayerSize,4)); # Contains biases for efficiency
global layerTwoWeights = rand!(zeros(hiddenLayerSize + 1)'); # Also contains biases

errorPlotData = open("cost.csv", "a");

for j = 1:typemax(UInt64) # Effectively infinite
        layerOneGradient = zeros(hiddenLayerSize,4);
        layerTwoGradient = zeros(hiddenLayerSize + 1)';

        error = 0;

        for i = 1:numOfDatapoints
                netOne = layerOneWeights*inputs[:,i];
                layerOneOutput = vcat(sigmoid.(netOne), 1); # Append 1 for bias term in weights

                netTwo = layerTwoWeights*layerOneOutput;
                layerTwoOutput = netTwo; # Identity output unit.

                cost = (layerTwoOutput - expectedOutput[i])^2;
                error += cost;
                dCost_dY = 2*(layerTwoOutput - expectedOutput[i]);

                layerOneGradient[:,:] += dCost_dY * layerTwoWeights[1:hiddenLayerSize] .* sigmoidDerivative.(netOne) .* inputs[:,i]';

                layerTwoGradient[:] += dCost_dY * layerOneOutput[:];
        end

        layerOneGradient /= numOfDatapoints;
        layerTwoGradient /= numOfDatapoints;

        global layerOneWeights -= learningRate * layerOneGradient;
        global layerTwoWeights -= learningRate * layerTwoGradient;

        print("Error: ", error/numOfDatapoints, "\n", "Num of iterations: ", j, '\n')

        if mod(j,50) == 1
                write(errorPlotData, string(j), ",", string(error/numOfDatapoints), "\n");
        end

        if mod(j,1000) == 0
                writedlm("layer1.csv", layerOneWeights, ',');
                writedlm("layer2.csv", layerTwoWeights, ',');
        end
end


text:
"""
    SampleSpace

Defines restrictions on the sample space of a two-stage design to be used during
optimization.
"""
mutable struct SampleSpace{TI<:Integer,TR<:Real}

  n1range::Vector{TI}
  nmax::TI
  n2min::TI
  maxnfact::TR
  nmincont::TI
  maxvariables::TI
  GS::Bool # group sequential ?
  specialnvalues::Vector{TI} # required for additional values always included in nvals
  specialcvalues::Vector{TI}

  function SampleSpace{TI,TR}(
    n1range::Vector{TI},
    nmax::TI,
    n2min::TI,
    maxnfact::TR,
    nmincont::TI,
    maxvariables::TI,
    GS::Bool
  ) where {TI<:Integer,TR<:Real}

    minimum(n1range) < 1    ? error("minimal n1 must be >= 1")    : nothing
    maximum(n1range) > nmax ? error("maximal n1 must be <= nmax") : nothing
    maxvariables < 100      ? error("maxvariables must be >= 100") : nothing
    maxnfact = (maxnfact == Inf) ? nmax / minimum(n1range) : maxnfact
    new(sort(n1range), nmax, n2min, maxnfact, nmincont, maxvariables, GS, convert(Vector{TI}, zeros(0)), convert(Vector{TI}, zeros(0)))

  end # inner constructor (SampleSpace)

end # SampleSpace

"""
    SampleSpace(
      n1range,
      nmax::TI;
      n2min::TI = 1,
      maxnfact::TR = Inf,
      nmincont::TI = 0,
      maxvariables::TI = 500000,
      GS::Bool = false
    ) where {TI<:Integer,TR<:Real}

Creates a [`SampleSpace`](@ref)-object which defines the search space for finding
optimal two-stage designs using [`optimaldesign`](@ref).
The primary use of dedicated sample sapce objects is the incoporation of
specific operational constraints on optimal designs (e.g. maximal sample size).

# Parameters

| Parameter      | Description |
| -------------: | :---------- |
| n1range        | possible stage-one sample sizes, can be anything which is convertable to an integer vector |
| nmax           | maximal overall sample size (stage one and two combined) |
| [n2min]        | minimal stage-two sample size |
| [maxnfact]     | maximal overall sample size must be smaller than `maxnfact * n1`  |
| [nmincont]     | minimal sample size upon rejection of the null hypothesis |
| [maxvariables] | (approximate) maximal number of binary variables used when finding optimal design; if the sample space is too big this can be used to find approximate solutions to the optimization problem |
| [GS]           | flag indicating wheather the design should be group-sequential (constant stage two sample size upon continuation) |
"""
function SampleSpace( # unspecific variant
  n1range,
  nmax::TI;
  n2min::TI = 1,
  maxnfact::TR = Inf,
  nmincont::TI = 0,
  maxvariables::TI = 500000,
  GS::Bool = false
) where {TI<:Integer,TR<:Real}

  SampleSpace(
    convert(Vector{TI}, n1range), nmax, n2min, maxnfact, nmincont, maxvariables, GS
  )

end # external constructor (SampleSpace)

function SampleSpace(
  n1range::Vector{TI},
  nmax::TI,
  n2min::TI,
  maxnfact::TR,
  nmincont::TI,
  maxvariables::TI,
  GS::Bool
) where {TI<:Integer,TR<:Real}

  SampleSpace{TI,TR}(
      n1range, nmax, n2min, maxnfact, nmincont, maxvariables, GS
  )

end # external constructor (SampleSpace)

Base.size(::SampleSpace) = ()

Base.getindex(ss::SampleSpace, i) = ss

Base.show(io::IO, ss::SampleSpace) = print("SampleSpace")

"""
    interimsamplesizerange(ss::SampleSpace)

Return integer vector of valid ``n_1`` values.
"""
interimsamplesizerange(ss::SampleSpace) = ss.n1range

"""
    maxsamplesize(ss::SampleSpace)

Return maximal overal sample size as integer.
"""
maxsamplesize(ss::SampleSpace) = ss.nmax

"""
    possible(
      n1::TI2, ss::SampleSpace{TI,TR}
    ) where {TI<:Integer,TR<:Real,TI2<:Integer}

Return Bool indicating wheather `n1` is possible under `ss`.
"""
function possible(
  n1::TI2, ss::SampleSpace{TI,TR}
) where {TI<:Integer,TR<:Real,TI2<:Integer}

  return n1 in ss.n1range

end # possible

"""
    possible(
      n1::TI, n::TI, c::TR, ss::SampleSpace{TI2,TR2}
    ) where {TI<:Integer,TR<:Real,TI2<:Integer,TR2<:Real}

Return Bool indicating wheather the combination of `n1`, `n`, and `c` is
possible under `ss`.
"""
function possible(
  n1::TI, n::TI, c::TR, ss::SampleSpace{TI2,TR2}
) where {TI<:Integer,TR<:Real,TI2<:Integer,TR2<:Real}

  res = n1 in ss.n1range # n1 must be possible
  res = n <= ss.nmax ? res : false # overall sample size constraint must not
                                     # be violated
  res = n1 <= n ? res : false # n1 cannot be larger than n
  res = n <= ss.maxnfact * n1 ? res : false # n cannot be too large in relation
                                            # to n1
  n2  = n - n1
  if (n2 < ss.n2min) & (c != Inf)
    # n2 too small (also affects "early stopping" for efficacy!)
    if (c == -Inf) & (n1 == n) & (n1 >= ss.nmincont) # fine, early stopping for efficacy
      nothing
    else
      res = false
    end
  end
  res = (n  < ss.nmincont) & (c != Inf) ? false : res # only stopping for futility may violate nmincount
  return res
end # possible

"""
    maxsamplesize(
      ss::SampleSpace{TI,TR}, n1::TI2
    ) where {TI<:Integer,TR<:Real,TI2<:Integer}

Return maximal sample size given `n1` as integer.
"""
function maxsamplesize(
  ss::SampleSpace{TI,TR}, n1::TI2
) where {TI<:Integer,TR<:Real,TI2<:Integer}

  possible(n1, ss)
  convert(TI, min(ss.nmax, floor(n1 * ss.maxnfact)))

end # maxsamplesize

"""
    isgroupsequential(ss::SampleSpace)

Indicate whether the sample space is restricted to group-sequential designs.
"""
isgroupsequential(ss::SampleSpace) = ss.GS

function getnvals(
  ss::SampleSpace{TI,TR}, n1::TI2
) where {TI<:Integer,TR<:Real,TI2<:Integer}

    nmax    = maxsamplesize(ss, n1)
    frac    = (nmax - n1 + 1) / sqrt(ss.maxvariables / n1)
    if frac > 1
        info(@sprintf("thinning nvals by factor %.3f", 1 / frac))
    else
        frac = 1
    end
    nvals   = convert(Vector{Int64}, round.(collect(n1:frac:nmax)))
    if !(n1 + ss.nmincont in nvals)
        push!(nvals, n1 + ss.nmincont)
    end
    if !((n1 + ss.n2min) in nvals)
        push!(nvals, n1 + ss.n2min)
    end
    if !(maxsamplesize(ss, n1) in nvals)
        push!(nvals, maxsamplesize(ss, n1))
    end
    for snval in ss.specialnvalues # ensure that all special nvals are contained!
        if !(snval in nvals) & (snval >= n1) & (snval <= nmax)
            push!(nvals, snval)
        end
    end
    return convert(Vector{Int64}, sort(nvals))

end # getnvals

function getcvals(
  ss::SampleSpace{TI,TR}, n1::TI2
) where {TI<:Integer,TR<:Real,TI2<:Integer}

    nmax        = maxsamplesize(ss, n1)
    frac        = nmax / sqrt(ss.maxvariables / n1)
    if frac > 1
        info(@sprintf("thinning cvals by factor %.3f", 1 / frac))
    else
        frac = 1
    end
    cvalsfinite = round.(collect(0:frac:(nmax - 1)))
    for scval in ss.specialcvalues # ensure that all special cvals are contained!
      if !(scval in cvalsfinite) & (scval >= n1) & (scval <= nmax)
          push!(cvalsfinite, scval)
      end
    end
    return cvalsfinite

end # getcvals

max_stars_repo_name: kkmann/binaryTwoStageDesigns

max_issues_repo_name: kkmann/binaryTwoStageDesigns

max_forks_repo_name: kkmann/binaryTwoStageDesigns


text:
import Oceananigans.Utils: prettytime
import Oceananigans.TimeSteppers: reset!

# It's not a model --- its a simulation!

default_progress(simulation) = nothing

mutable struct Simulation{ML, TS, DT, ST, DI, OW, CB}
              model :: ML
        timestepper :: TS
                 Δt :: DT
     stop_iteration :: Float64
          stop_time :: ST
    wall_time_limit :: Float64
        diagnostics :: DI
     output_writers :: OW
          callbacks :: CB
      run_wall_time :: Float64
            running :: Bool
        initialized :: Bool
end

"""
    Simulation(model; Δt,
               stop_iteration = Inf,
               stop_time = Inf,
               wall_time_limit = Inf)

Construct a `Simulation` for a `model` with time step `Δt`.

Keyword arguments
=================

- `Δt`: Required keyword argument specifying the simulation time step. Can be a `Number`
        for constant time steps or a `TimeStepWizard` for adaptive time-stepping.

- `stop_iteration`: Stop the simulation after this many iterations.

- `stop_time`: Stop the simulation once this much model clock time has passed.

- `wall_time_limit`: Stop the simulation if it's been running for longer than this many
                     seconds of wall clock time.
"""
function Simulation(model; Δt,
                    stop_iteration = Inf,
                    stop_time = Inf,
                    wall_time_limit = Inf)

   if stop_iteration == Inf && stop_time == Inf && wall_time_limit == Inf
       @warn "This simulation will run forever as stop iteration = stop time " *
             "= wall time limit = Inf."
   end

   diagnostics = OrderedDict{Symbol, AbstractDiagnostic}()
   output_writers = OrderedDict{Symbol, AbstractOutputWriter}()
   callbacks = OrderedDict{Symbol, Callback}()

   callbacks[:stop_time_exceeded] = Callback(stop_time_exceeded)
   callbacks[:stop_iteration_exceeded] = Callback(stop_iteration_exceeded)
   callbacks[:wall_time_limit_exceeded] = Callback(wall_time_limit_exceeded)

   # Check for NaNs in the model's first prognostic field every 100 iterations.
   model_fields = fields(model)
   field_to_check_nans = NamedTuple{keys(model_fields) |> first |> tuple}(first(model_fields) |> tuple)
   nan_checker = NaNChecker(field_to_check_nans)
   callbacks[:nan_checker] = Callback(nan_checker, IterationInterval(100))

   # Convert numbers to floating point; otherwise preserve type (eg for DateTime types)
   FT = eltype(model.grid)
   Δt = Δt isa Number ? FT(Δt) : Δt
   stop_time = stop_time isa Number ? FT(stop_time) : stop_time

   return Simulation(model,
                     model.timestepper,
                     Δt,
                     Float64(stop_iteration),
                     stop_time,
                     Float64(wall_time_limit),
                     diagnostics,
                     output_writers,
                     callbacks,
                     0.0,
                     false,
                     false)
end

Base.show(io::IO, s::Simulation) =
    print(io, "Simulation{$(typeof(s.model).name){$(Base.typename(typeof(s.model.architecture))), $(eltype(s.model.grid))}}\n",
              "├── Model clock: time = $(prettytime(s.model.clock.time)), iteration = $(s.model.clock.iteration)", '\n',
              "├── Next time step: $(prettytime(s.Δt))", '\n',
              "├── Elapsed wall time: $(prettytime(s.run_wall_time))", '\n',
              "├── Stop time: $(prettytime(s.stop_time))", '\n',
              "├── Stop iteration : $(s.stop_iteration)", '\n',
              "├── Wall time limit: $(s.wall_time_limit)", '\n',
              "├── Callbacks: $(ordered_dict_show(s.callbacks, "│"))", '\n',
              "├── Output writers: $(ordered_dict_show(s.output_writers, "│"))", '\n',
              "└── Diagnostics: $(ordered_dict_show(s.diagnostics, "│"))")

#####
##### Utilities
#####

"""
    time(sim::Simulation)

Return the current simulation time.
"""
Base.time(sim::Simulation) = sim.model.clock.time

"""
    iteration(sim::Simulation)

Return the current simulation iteration.
"""
iteration(sim::Simulation) = sim.model.clock.iteration

"""
    prettytime(sim::Simulation)

Return `sim.model.clock.time` as a prettily formatted string."
"""
prettytime(sim::Simulation, longform=true) = prettytime(time(sim))

"""
    run_wall_time(sim::Simulation)

Return `sim.run_wall_time` as a prettily formatted string."
"""
run_wall_time(sim::Simulation) = prettytime(sim.run_wall_time)

"""
    reset!(sim)

Reset `sim`ulation, `model.clock`, and `model.timestepper` to their initial state.
"""
function reset!(sim::Simulation)
    sim.model.clock.time = 0.0
    sim.model.clock.iteration = 0
    sim.model.clock.stage = 1
    sim.stop_iteration = Inf
    sim.stop_time = Inf
    sim.wall_time_limit = Inf
    sim.initialized = false
    sim.running = true
    reset!(sim.model.timestepper)
    return nothing
end

#####
##### Default stop criteria callback functions
#####

function stop_iteration_exceeded(sim)
    if sim.model.clock.iteration >= sim.stop_iteration
        @info "Simulation is stopping. Model iteration $(sim.model.clock.iteration) " *
               "has hit or exceeded simulation stop iteration $(Int(sim.stop_iteration))."
       sim.running = false
    end
    return nothing
end

function stop_time_exceeded(sim)
    if sim.model.clock.time >= sim.stop_time
       @info "Simulation is stopping. Model time $(prettytime(sim.model.clock.time)) " *
             "has hit or exceeded simulation stop time $(prettytime(sim.stop_time))."
       sim.running = false
    end
    return nothing
end

function wall_time_limit_exceeded(sim)
    if sim.run_wall_time >= sim.wall_time_limit
        @info "Simulation is stopping. Simulation run time $(run_wall_time(sim)) " *
              "has hit or exceeded simulation wall time limit $(prettytime(sim.wall_time_limit))."
       sim.running = false
    end
    return nothing
end


text:
include("empty.jl")
include("fourrooms.jl")
include("gotodoor.jl")


text:
"""
    make_data(n_data::Int, p::Vector{Float64})

Generate some example dyadic data using a Gaussian mixture model.
"""
function make_data(n_data::Int, p::Vector{Float64})

    @assert n_data >= 2
    @assert sum(p) == 1.0
    @assert all(p .>= 0)
    @assert length(p) == 3

    c = cumsum(p)
    r = rand(n_data)
    data_A = -1.0 * (r .<= c[1]) + 1.0 * (r .>= c[2])
    data_V = randn((n_data, n_data))
    data = UpperTriangular(data_A .* data_A' .+ data_V)

    return data
end

"""
    mean(x::Vector{<:Number})

Compute the mean of a vector of numbers.
"""
function mean(x::Vector{<:Number})

    return sum(x) / length(x)
end

"""
    phi(t::Real)

Compute the standard normal density function.
"""
function phi(t::Real)

    return (2 * pi)^(-0.5) * exp(-(t^2) / 2)
end

"""
    get_f(est::DyadicKernelDensityEstimator)

Get the true density function from example dyadic Gaussian mixture data.
"""
function get_f(est::DyadicKernelDensityEstimator)

    p = est.meta["p"]
    f = (p[1]^2 + p[3]^2) * phi.(est.evals .- 1)
    f += p[2]*(2 - p[2]) * phi.(est.evals)
    f += 2*p[1]*p[3] * phi.(est.evals .+ 1)
    return f
end

"""
    get_RIMSE(est::DyadicKernelDensityEstimator)

Compute the root integrated mean squared error.
"""
function get_RIMSE(est::DyadicKernelDensityEstimator)

    f = get_f(est)
    RMSE = sqrt(mean((est.fhat .- f).^2))
    return RMSE
end

"""
    get_ucb_coverage(est::DyadicKernelDensityEstimator)

Check if the uniform confidence band covers the true density function.
"""
function get_ucb_coverage(est::DyadicKernelDensityEstimator)

    return all(est.ucb[1,:] .<= get_f(est) .<= est.ucb[2,:])
end

"""
    get_pci_coverage(est::DyadicKernelDensityEstimator)

Check if the pointwise confidence intervals all cover the true density function.
"""
function get_pci_coverage(est::DyadicKernelDensityEstimator)

    return all(est.pci[1,:] .<= get_f(est) .<= est.pci[2,:])
end

"""
    get_bci_coverage(est::DyadicKernelDensityEstimator)

Check if the Bonferroni confidence intervals all cover the true density function.
"""
function get_bci_coverage(est::DyadicKernelDensityEstimator)

    return all(est.bci[1,:] .<= get_f(est) .<= est.bci[2,:])
end

"""
    get_ucb_average_width(est::DyadicKernelDensityEstimator)

Return the average width of the uniform confidence band.
"""
function get_ucb_average_width(est::DyadicKernelDensityEstimator)

    return mean(est.ucb[2,:] .- est.ucb[1,:])
end

"""
    get_pci_average_width(est::DyadicKernelDensityEstimator)

Return the average width of the pointwise confidence intervals.
"""
function get_pci_average_width(est::DyadicKernelDensityEstimator)

    return mean(est.pci[2,:] .- est.pci[1,:])
end

"""
    get_bci_average_width(est::DyadicKernelDensityEstimator)

Return the average width of the Bonferroni confidence intervals.
"""
function get_bci_average_width(est::DyadicKernelDensityEstimator)

    return mean(est.bci[2,:] .- est.bci[1,:])
end


text:
#
# This file is a part of MolecularGraph.jl
# Licensed under the MIT License http://opensource.org/licenses/MIT
#

module MolecularGraphUtilTest
    using Test
    using MolecularGraph.Util
    include("./util/iterator.jl")
    include("./util/math.jl")
end

module MolecularGraphGraphTest
    using Test
    using MolecularGraph.Graph

    include("./graph/dfs.jl")
    include("./graph/shortestpath.jl")
    include("./graph/generator.jl")

    include("./graph/plaingraph.jl")
    include("./graph/plaindigraph.jl")
    include("./graph/plainhypergraph.jl")

    include("./graph/multigraph.jl")
    include("./graph/dag.jl")
    include("./graph/cycle.jl")
    include("./graph/connectivity.jl")
    include("./graph/inducedsubgraph.jl")
    include("./graph/disjointunion.jl")
    include("./graph/linegraph.jl")

    include("./graph/triangle.jl")
    include("./graph/clique.jl")
    include("./graph/bipartite.jl")
    include("./graph/planarity.jl")
    include("./graph/product.jl")
    include("./graph/isomorphism/vf2.jl")
    include("./graph/isomorphism/cliquebased.jl")
end

module MolecularGraphGeometryTest
    using Test
    using LinearAlgebra
    using MolecularGraph.Geometry

    include("./geometry/cartesian.jl")
    include("./geometry/internal.jl")
end

module MolecularGraphTest
    using Test
    using MolecularGraph
    using MolecularGraph.Graph

    include("./model/atom.jl")
    include("./model/molgraph.jl")

    include("sdfilereader.jl")
    include("sdfilewriter.jl")

    include("./smarts/base.jl")
    include("./smarts/logicaloperator.jl")
    include("./smarts/atom.jl")
    include("./smarts/bond.jl")
    include("./smarts/molecule.jl")
    include("./smarts/smiles.jl")
    include("./smarts/smarts.jl")

    include("properties.jl")
    include("preprocess.jl")
    include("stereo.jl")
    include("mass.jl")
    include("wclogp.jl")
    include("substructure.jl")
    include("funcgroup.jl")
    include("mcs.jl")
    include("inchi.jl")

    include("./draw/base.jl")
    include("./draw/svg.jl")
end


text:
# A collection of functions to load and preprocess data for training.
module DataHelpers
using MLDatasets
using Base.Iterators: partition
using Statistics
using Random

################################################################################
# Dataset partitioning
################################################################################
function makebatches(data, labels; batchsize=128)
    batchinds = collect(partition(1:length(labels), batchsize))
    n = ndims(data)
    databatches = map(batch -> data[repeat([Colon()], n-1)..., batch], batchinds)
    labelbatches = map(batch -> labels[batch], batchinds)
    return databatches, labelbatches
end # function makebatches

################################################################################
# Preprocessing
################################################################################
"""
    standardize(imgs)

Normalize each image in `imgs` to have mean 0 and variance 1. Assumes that the
first two dimensions are the image dimensions.

!!! warning
    This should only be used on grayscale images. Color channels would be treated
    like independent images.
"""
function standardize(imgs)
    imgs = imgs .- mean(imgs, dims=(1,2))
    imgs = imgs ./ std(imgs, dims=(1,2))
    return imgs
end # function normalize_images

"""
    whitenoise(imgs, σ)

Add white noise with variance `σ` to `imgs`.
"""
whitenoise(imgs, σ) = return imgs + randn(eltype(imgs), size(imgs)...) .* σ

"""
    saltpepper(imgs, p)

Add salt-and-pepper noise to 'imgs' by randomly swithing off pixels with probability
`p` (half will be set to 0, half to 1).
"""
function saltpepper(imgs, p)
    imgs = copy(imgs)
    inds = shuffle(randsubseq(1:length(imgs), 0.1))
    imgs[inds[1:2:end]] .= 0.0
    imgs[inds[2:2:end]] .= 1.0
    return imgs
end

################################################################################
# Dataset generators
################################################################################
function makemultiMNIST(; set=:training, digits=2, offset=4)
    if set == :training
        imgs, lbls = MNIST.traindata()
    elseif set == :test
        imgs, lbls = MNIST.testdata()
    end
    imgs = Float64.(imgs)
    x = size(imgs, 1)
    y = size(imgs, 2)
    n = size(imgs, 3)
    newimgs = zeros(x+2offset, y+2offset, n)
    inds = [randperm(n) for i in 1:digits]
    for i in 1:n
        for j in 1:digits
            Δx, Δy = rand(-offset:offset, 2)
            newimgs[1+offset-Δx:x+offset-Δx, 1+offset-Δy:y+offset-Δy, i] = max.(
                imgs[:, :, inds[j][i]],
                newimgs[1+offset-Δx:x+offset-Δx, 1+offset-Δy:y+offset-Δy, i]
            )
        end
    end
    newlbls = [lbls[inds[i]] for i in 1:digits]
    newlbls = collect(zip(newlbls...))
    return newimgs, newlbls
end # function makemultiMNIST
end # module DataHelpers


text:
module Tanks

export M4Sherman

"M4Sherman is a pivotal WW2 tank from the US."
mutable struct M4Sherman

    "power status: true = on, false = off"
    power::Bool

    "current direction in radians"
    direction::Float64

    "current position coordinate (x,y)"
    position::Tuple{Float64, Float64}

    wheels::Tuple{Symbol, Symbol}

end

# Import generic functions
import Vehicle: power_on!, power_off!, turn!, move!, position, engage_wheels!, has_wheels

# Implementation of Vehicle interface

function power_on!(t::M4Sherman)
    t.power = true
    println("Powered on: ", t, " ready to go!")
    nothing
end

function power_off!(t::M4Sherman)
    t.power = false
    println("Powered off: ", t)
    nothing
end

function turn!(t::M4Sherman, direction)
    t.direction = direction
    println("Now barreling towards ", direction, ": ", t)
    nothing
end

function move!(t::M4Sherman, distance)
    x, y = t.position
    dx = round(distance * cos(t.direction), digits = 2)
    dy = round(distance * sin(t.direction), digits = 2)
    t.position = (x + dx, y + dy)
    println("Moved (", dx, ",", dy, "): ", t)
    nothing
end

function position(t::M4Sherman)
    t.position
end

has_wheels(t::M4Sherman) = true

function engage_wheels!(t::M4Sherman)
    println("Making mark with both wheels")
    nothing
end

end # module


text:
using ReachabilityBenchmarks, MathematicalSystems, Reachability, Plots

include(@relpath "beam_model.jl")
include(@relpath "beam_specifications.jl")

S = beam_model()
X0, options = beam_specification()
options = Dict(options)

# initial value problem
problem = InitialValueProblem(S, X0)

# partition with 1D blocks
partition = [i:i for i in 1:348]

# reachability algorithms
algorithm_reach = BFFPSV18(:δ => 1e-3, :vars => [89], :partition => partition)
algorithm_check = BFFPSV18(:δ => 5e-5, :vars => [89], :partition => partition)

# check property
options[:mode] = "check"
solution = solve(problem, options; op=algorithm_check)
@assert solution.satisfied

# compute flowpipe
options[:mode] = "reach"
solution = solve(problem, options; op=algorithm_reach)

# project flowpipe to time and x₈₉
solution.options[:plot_vars] = [0, 89]
solution_proj = project(solution)

# plot projection
plot(solution_proj)
savefig("beam")

max_stars_repo_name: JuliaReach/ReachabilityBenchmarks

max_issues_repo_name: JuliaReach/ReachabilityBenchmarks

max_forks_repo_name: JuliaReach/ReachabilityBenchmarks


text:
"""
    save(f::FileIO.File)
    save(f::FileIO.File, data, lidict)
    save(filename::AbstractString, data, ldict)

Save profiling data to a file. If `data` and `lidict` are not supplied, they are obtained
from

    data, lidict = Profile.retrieve()

Note that the data saved to the file discard some system-specific information to allow portability.
Some visualization modes, like [`StackFrameCategory`](@ref), are not available for data loaded from
such files.

These files conventionally have the extension ".jlprof".
If you just supply a string filename ending with this extension,
you must pass `data` and `lidict` explicitly, because
FileIO has its own interpretation of the meaning of `save` with no arguments.
"""
function save(f::File{format"JLPROF"}, data::AbstractVector{<:Unsigned}, lidict::Profile.LineInfoDict)
    data_u64 = convert(AbstractVector{UInt64}, data)
    open(f, "w") do io
        write(io, magic(format"JLPROF"))
        # Write an endianness revealer
        write(io, 0x01020304)
        # Write an indicator that this is data/lidict format
        write(io, 0x01)
        write(io, Int64(length(data_u64)))
        write(io, data_u64)
        write(io, Int64(length(lidict)))
        for (k, v) in lidict
            write(io, k)
            write(io, Int32(length(v)))
            for sf in v
                symwrite(io, sf.func)
                symwrite(io, sf.file)
                write(io, sf.line)
                write(io, sf.from_c)
                write(io, sf.inlined)
            end
        end
    end
    return nothing
end

# Note: this doesn't work from FileIO because it returns an anonymous function for saving data
save(f::File{format"JLPROF"}) = save(f, Profile.retrieve()...)

"""
    data, lidict = load(f::FileIO.File)
    data, lidict = load(filename::AbstractString)

Load profiling data. You can reconstruct the flame graph from `flamegraph(data; lidict=lidict)`.
"""
function load(f::File{format"JLPROF"})
    open(f) do io
        skipmagic(io)
        load(io)
    end
end

function load(io::Stream{format"JLPROF"})
    endian = read(io, UInt32)
    endian == 0x01020304 || error("bswap not yet supported, please report as an issue to FlameGraphs.jl")
    fmt = read(io, UInt8)
    if fmt == 0x01
        # This is data/lidict format
        n = read(io, Int64)
        data = Vector{UInt64}(undef, n)
        read!(io, data)
        n = read(io, Int64)
        lidict = Profile.LineInfoDict()
        for i = 1:n
            k = read(io, UInt64)
            nsf = read(io, Int32)
            sfs = StackFrame[]
            for j = 1:nsf
                func    = symread(io)
                file    = symread(io)
                line    = read(io, Int)
                from_c  = read(io, Bool)
                inlined = read(io, Bool)
                push!(sfs, StackFrame(func, file, line, nothing, from_c, inlined, 0x0))
            end
            lidict[k] = sfs
        end
        return data, lidict
    end
    error("format ", fmt, " not recognized")
end

function symwrite(io, sym::Symbol)
    str = string(sym)
    write(io, Int32(ncodeunits(str)))
    write(io, str)
end

function symread(io)
    n = read(io, Int32)
    b = read(io, n)
    return Symbol(String(b))
end


text:
@testset "ultrametrize!" begin
n = readTopology("((((B)#H1:::0.7,D):4.0,(#H1:::0.3,E):6.2):2.0,O);")
@test @test_logs QuartetNetworkGoodnessFit.ultrametrize!(n, true)
@test writeTopology(n, round=true) == "((((B:0.0)#H1:2.2::0.7,D:2.2):4.0,(#H1:0.0::0.3,E:0.0):6.2):2.0,O:8.2);"
n = readTopology("((((B)#H1:::0.7,D):4.0,(#H1:1.0::0.3,E):2.2):2.0,O);")
@test !(@test_logs (:warn, r"^could not adjust length of major edge") QuartetNetworkGoodnessFit.ultrametrize!(n, true))
@test writeTopology(n, round=true) == "((((B:0.0)#H1:0.0::0.7,D:0.0):4.0,(#H1:1.0::0.3,E:1.8):2.2):2.0,O:6.0);"
n = readTopology("((((B:2.0)#H1:0.5::0.7,D):4.0,(#H1:::0.3,E):6.2):2.0,O);")
@test !(@test_logs (:warn, r"^could not adjust length of minor edge") QuartetNetworkGoodnessFit.ultrametrize!(n, true))
@test writeTopology(n, round=true) == "((((B:2.0)#H1:0.5::0.7,D:2.5):4.0,(#H1:0.0::0.3,E:0.3):6.2):2.0,O:8.5);"
n = readTopology("((((B)#H1:0.6::0.7,D):4.0,(#H1:0.1::0.3,E):6.2):2.0,O);")
@test !(@test_logs (:warn,r"hybrid H") QuartetNetworkGoodnessFit.ultrametrize!(n, true))
end

@testset "reroot!" begin
rn = readTopology("((((D:0.4,C:0.4):4.8,((A:0.8,B:0.8):2.2)#H1:2.2::0.7):4.0,(#H1:0::0.3,E:3.0):6.2):2.0,O:11.2);");
tn = readTopology("((D:0.4,C:0.4):4.8,((A:0.8,B:0.8):2.2)#H1:2.2::0.7,((#H1:0.0::0.3,E:3.0):6.2,(O:11.2):2.0):4.0);")
@test QuartetNetworkGoodnessFit.reroot!(tn, rn) == 0
tn = readTopology("((D:0.4,C:0.4):4.8,((A:0.8,B:0.8):2.2)#H1:2.2::0.7,((#H1:0.0::0.3,E:3.0):6.2,O:13.2):4.0);")
# same as before, except for degree-2 node removed
@test QuartetNetworkGoodnessFit.reroot!(tn, rn) == 1
end


text:
function part1()
    lines = readlines("inputD2.txt")
    horiz = 0
    depth = 0
    for line in lines
        instruction, amount = split(line, " ")
        amount = parse(Int, amount)

        if instruction == "forward"
            horiz += amount
        elseif instruction == "up"
            depth -= amount
        elseif instruction == "down"
            depth += amount
        end
    end

    return horiz * depth
end

function part2()
    lines = readlines("inputD2.txt")
    horiz = 0
    depth = 0
    aim = 0
    for line in lines
        instruction, amount = split(line, " ")
        amount = parse(Int, amount)

        if instruction == "forward"
            horiz += amount
            depth += aim * amount
        elseif instruction == "up"
            aim -= amount
        elseif instruction == "down"
            aim += amount
        end
    end

    return horiz * depth
end

println(part1())
println(part2())


text:
@testset "FlatZinc" begin
    include("model.jl")
    include("export.jl")
    include("import.jl")
    include("optimizer.jl")
end


text:
# automatically generated by the FlatBuffers compiler, do not modify

include("UnionVector/Attacker.jl")
include("UnionVector/Rapunzel.jl")
include("UnionVector/BookReader.jl")
include("UnionVector/Character.jl")
include("UnionVector/Movie.jl")

max_stars_repo_name: rjkat/flatbuffers-julia

max_issues_repo_name: rjkat/flatbuffers-julia

max_forks_repo_name: rjkat/flatbuffers-julia


text:
module Met_Basic_Types

function RemoveTr(xarr :: Array{Float64,2})
  time = ones(2,size(xarr,2));
  time[2,:] = 0.0 + [1:size(xarr,2);];
  xarr = xarr - (xarr/time)*time;
  xarr = broadcast(-,xarr, mean(xarr,2))
end
export RemoveTr

# package code goes here
include("Areas/Areas.jl")
include("Model_Operations/area.jl")
include("Times/Time_Operator.jl")
include("Time_1D_Operators/All_Timers.jl")
include("Dimen_Forward_and_Backward/All_Forwardors.jl")
include("often_used/all.jl")
end # module


text:
push!(LOAD_PATH, "/home/ian/Grad/Research/Disks/DiskJockey/")

using ArgParse

s = ArgParseSettings()
@add_arg_table s begin
    # "--opt1"
    # help = "an option with an argument"
    # default = 0
    "--norad"
    help = "Use the image already here."
    action = :store_true
    "config"
    help = "a YAML configuration file"
    required = true
end

parsed_args = parse_args(ARGS, s)

import YAML
config = YAML.load(open(parsed_args["config"]))

# Make some diagnostic plots to show what the model looks like in analytic terms

using emodel
import PyPlot.plt
using LaTeXStrings
using constants

pp = config["parameters"]
params = ["M_star", "a_c", "T_10", "q", "gamma", "logM_CO", "ksi", "dpc", "incl", "PA", "e", "w", "vel", "mu_RA", "mu_DEC"]
nparam = length(params)
starting_param = Array(Float64, nparam)

for i=1:nparam
    starting_param[i] = pp[params[i]][1]
end

# Convert logM_CO to M_CO
starting_param[6] = 10^starting_param[6]
pars = Parameters(starting_param...)

grd = config["grid"]
grid = Grid(grd["nr"], grd["ntheta"], grd["nphi"], grd["r_in"], grd["r_out"], grd["na"], true)

# Plot the rings and ring centers for the disk.
function plot_rings(pars::Parameters, grid)
    fig = plt.figure(figsize=(8,8))
    ax = fig[:add_subplot](111, polar=true)

    phis = linspace(0, 2pi, 100)

    # Plot the grid cell edges in radius and phi
    for R in grid.Rs
        rr = R * ones(phis)
        ax[:plot](phis, rr ./AU, "b", lw=0.1)
    end

    for Phi in grid.Phis
        ax[:plot](Phi * ones(grid.Rs), grid.Rs ./AU, "b", lw=0.1)
    end

    # theta, r
    for a in grid.as
        rr = emodel.radius(phis, a, pars)
        ax[:plot](phis, rr./AU, "k-.")
    end

    for A in grid.As
        RR = emodel.radius(phis, A, pars)
        ax[:plot](phis, RR./AU, "k-")
    end

    # Plot the outermost wall radii
    rout = grid.Rs[end] .* ones(phis)
    ax[:plot](phis, rout./AU, "b-")

    fig[:subplots_adjust](left=0.15, bottom=0.15, right=0.85)

    plt.savefig("../plots/radii_outer.png")
end

function plot_inner_rings(pars::Parameters, grid)
    fig = plt.figure(figsize=(8,8))
    ax = fig[:add_subplot](111, polar=true)

    phis = linspace(0, 2pi, 100)

    # Plot the grid cell edges in radius and phi
    for R in grid.Rs
        rr = R * ones(phis)
        ax[:plot](phis, rr ./AU, "b", lw=0.1)
    end

    for Phi in grid.Phis
        ax[:plot](Phi * ones(grid.Rs), grid.Rs ./AU, "b", lw=0.1)
    end

    # theta, r
    for a in grid.as
        rr = emodel.radius(phis, a, pars)
        ax[:plot](phis, rr./AU, "k-.")
    end

    for A in grid.As
        RR = emodel.radius(phis, A, pars)
        ax[:plot](phis, RR./AU, "k-")
    end

    # Plot the innermost radii
    rin = grid.Rs[1] .* ones(phis)
    ax[:plot](phis, rin./AU, "b-")

    ax[:set_rmax](2)

    fig[:subplots_adjust](left=0.15, bottom=0.15, right=0.85)

    plt.savefig("../plots/radii_inner.png")
end

plot_rings(pars, grid)
plot_inner_rings(pars, grid)

max_stars_repo_name: elnjensen/DiskJockey

max_issues_repo_name: elnjensen/DiskJockey

max_forks_repo_name: elnjensen/DiskJockey


text:
@testset "mh.jl" begin
    @turing_testset "mh constructor" begin
        Random.seed!(0)
        N = 500
        s1 = MH(
            (:s, InverseGamma(2,3)),
            (:m, GKernel(3.0)))
        s2 = MH(:s, :m)
        s3 = MH()
        for s in (s1, s2, s3)
            @test DynamicPPL.alg_str(Sampler(s, gdemo_default)) == "MH"
        end

        c1 = sample(gdemo_default, s1, N)
        c2 = sample(gdemo_default, s2, N)
        c3 = sample(gdemo_default, s3, N)

        s4 = Gibbs(MH(:m), MH(:s))
        c4 = sample(gdemo_default, s4, N)
    end
    @numerical_testset "mh inference" begin
        Random.seed!(125)
        alg = MH()
        chain = sample(gdemo_default, alg, 2000)
        check_gdemo(chain, atol = 0.1)

        # MH with Gaussian proposal
        alg = MH(
            (:s, InverseGamma(2,3)),
            (:m, GKernel(1.0)))
        chain = sample(gdemo_default, alg, 7000)
        check_gdemo(chain, atol = 0.1)

        # MH within Gibbs
        alg = Gibbs(MH(:m), MH(:s))
        chain = sample(gdemo_default, alg, 2000)
        check_gdemo(chain, atol = 0.1)

        # MoGtest
        gibbs = Gibbs(
            CSMC(15, :z1, :z2, :z3, :z4),
            MH((:mu1,GKernel(1)), (:mu2,GKernel(1)))
        )
        chain = sample(MoGtest_default, gibbs, 5000)
        check_MoGtest_default(chain, atol = 0.15)
    end

    # Test MH shape passing.
    @turing_testset "shape" begin
        @model M(mu, sigma, observable) = begin
            z ~ MvNormal(mu, sigma)

            m = Array{Float64}(undef, 1, 2)
            m[1] ~ Normal(0, 1)
            m[2] ~ InverseGamma(2, 1)
            s ~ InverseGamma(2, 1)

            observable ~ Bernoulli(cdf(Normal(), z' * z))

            1.5 ~ Normal(m[1], m[2])
            -1.5 ~ Normal(m[1], m[2])

            1.5 ~ Normal(m[1], s)
            2.0 ~ Normal(m[1], s)
        end

        model = M(zeros(2), ones(2), 1)
        sampler = Inference.Sampler(MH(), model)

        dt, vt = Inference.dist_val_tuple(sampler, Turing.VarInfo(model))

        @test dt[:z] isa AdvancedMH.StaticProposal{<:MvNormal}
        @test dt[:m] isa AdvancedMH.StaticProposal{Vector{ContinuousUnivariateDistribution}}
        @test dt[:m].proposal[1] isa Normal && dt[:m].proposal[2] isa InverseGamma
        @test dt[:s] isa AdvancedMH.StaticProposal{<:InverseGamma}

        @test vt[:z] isa Vector{Float64} && length(vt[:z]) == 2
        @test vt[:m] isa Vector{Float64} && length(vt[:m]) == 2
        @test vt[:s] isa Float64

        chain = sample(model, MH(), 100)

        @test chain isa MCMCChains.Chains
    end

    @turing_testset "proposal matrix" begin
        Random.seed!(100)

        mat = [1.0 -0.05; -0.05 1.0]

        prop1 = mat # Matrix only constructor
        prop2 = AdvancedMH.RandomWalkProposal(MvNormal(mat)) # Explicit proposal constructor

        spl1 = MH(prop1)
        spl2 = MH(prop2)

        # Test that the two constructors are equivalent.
        @test spl1.proposals.proposal.μ == spl2.proposals.proposal.μ
        @test spl1.proposals.proposal.Σ.mat == spl2.proposals.proposal.Σ.mat

        # Test inference.
        chain1 = sample(gdemo_default, spl1, 10000)
        chain2 = sample(gdemo_default, spl2, 10000)

        check_gdemo(chain1)
        check_gdemo(chain2)
    end

    @turing_testset "gibbs MH proposal matrix" begin
        # https://github.com/TuringLang/Turing.jl/issues/1556

        # generate data
        x = rand(Normal(5, 10), 20)
        y = rand(LogNormal(-3, 2), 20)

        # Turing model
        @model function twomeans(x, y)
            # Set Priors
            μ ~ MvNormal(2, 3)
            σ ~ filldist(Exponential(1), 2)

            # Distributions of supplied data
            x .~ Normal(μ[1], σ[1])
            y .~ LogNormal(μ[2], σ[2])

        end
        mod = twomeans(x, y)

        # generate covariance matrix for RWMH
        # with small-valued VC matrix to check if we only see very small steps
        vc_μ = convert(Array, 1e-4*I(2))
        vc_σ = convert(Array, 1e-4*I(2))

        alg = Gibbs(
            MH((:μ, vc_μ)),
            MH((:σ, vc_σ)),
        )

        chn = sample(
            mod,
            alg,
            3_000 # draws
        )

        chn2 = sample(mod, MH(), 3_000)

        # Test that the small variance version is actually smaller.
        v1 = var(diff(Array(chn["μ[1]"]), dims=1))
        v2 = var(diff(Array(chn2["μ[1]"]), dims=1))

        # FIXME: Do this properly. It sometimes fails.
        # @test v1 < v2
    end

    @turing_testset "vector of multivariate distributions" begin
        @model function test(k)
            T = Vector{Vector{Float64}}(undef, k)
            for i in 1:k
                T[i] ~ Dirichlet(5, 1.0)
            end
        end

        Random.seed!(100)
        chain = sample(test(1), MH(), 5_000)
        for i in 1:5
            @test mean(chain, "T[1][$i]") ≈ 0.2 atol=0.01
        end

        Random.seed!(100)
        chain = sample(test(10), MH(), 5_000)
        for j in 1:10, i in 1:5
            @test mean(chain, "T[$j][$i]") ≈ 0.2 atol=0.01
        end
    end

    @turing_testset "MH link/invlink" begin
        vi_base = DynamicPPL.VarInfo(gdemo_default)

        # Don't link when no proposals are given since we're using priors
        # as proposals.
        vi = deepcopy(vi_base)
        alg = MH()
        spl = DynamicPPL.Sampler(alg)
        Turing.Inference.maybe_link!(vi, spl, alg.proposals)
        @test !DynamicPPL.islinked(vi, spl)

        # Link if proposal is `AdvancedHM.RandomWalkProposal`
        vi = deepcopy(vi_base)
        d = length(vi_base[DynamicPPL.SampleFromPrior()])
        alg = MH(AdvancedMH.RandomWalkProposal(MvNormal(d, 1.0)))
        spl = DynamicPPL.Sampler(alg)
        Turing.Inference.maybe_link!(vi, spl, alg.proposals)
        @test DynamicPPL.islinked(vi, spl)

        # Link if ALL proposals are `AdvancedHM.RandomWalkProposal`.
        vi = deepcopy(vi_base)
        alg = MH(:s => AdvancedMH.RandomWalkProposal(Normal()))
        spl = DynamicPPL.Sampler(alg)
        Turing.Inference.maybe_link!(vi, spl, alg.proposals)
        @test DynamicPPL.islinked(vi, spl)

        # Don't link if at least one proposal is NOT `RandomWalkProposal`.
        # TODO: make it so that only those that are using `RandomWalkProposal`
        # are linked! I.e. resolve https://github.com/TuringLang/Turing.jl/issues/1583.
        # https://github.com/TuringLang/Turing.jl/pull/1582#issuecomment-817148192
        vi = deepcopy(vi_base)
        alg = MH(
            :m => AdvancedMH.StaticProposal(Normal()),
            :s => AdvancedMH.RandomWalkProposal(Normal())
        )
        spl = DynamicPPL.Sampler(alg)
        Turing.Inference.maybe_link!(vi, spl, alg.proposals)
        @test !DynamicPPL.islinked(vi, spl)
    end
end


text:
using Revise, Pkg

if  pwd() ∉ LOAD_PATH
  push!(LOAD_PATH,pwd())
end

Pkg.activate(pwd())
using LeverageChanges

#TODO: #1) Try scaling the returns by 100
#2) dup the first few rows of Hartzmark pg 652

#LeverageChanges.recompressfiles()
#LeverageChanges.recompressfile("CRSP-D")
#LeverageChanges.recompressfile("CRSP-M")

#refreshcomp (PreprocessCOMP): refreshes the compustat data from the csv.lz4 file
#refreshcrsp (PreprocessCRSP): refreshes the crsp data from the csv.lz4 file
#refreshmerge[1] (MergeCRSPCOMP): refreshes the compustat-crsp merge.
#refreshfactors (FormFactors): refreshes the compustat-crsp merge.
#[1]True if refreshcomp, refreshcrsp, or validatemerged are true
#[2]True if refreshcomp, refreshcrsp, refreshmerge, or validatemerged are true
#hello
#usequarterlydates
#

function leveragechangesscript()
  #LeverageChanges.recompressfile("CCM-A")
  #LeverageChanges.recompressfile("CRSP-D-I")
  #LeverageChanges.recompressfile("CRSP-D-II")
  #LeverageChanges.recompressfile("CRSP-M")
  #LeverageChanges.recompressfile("CCM")
  #LeverageChanges.recompressfile("COMP-A")

  SOUP_2_NUTS::Bool = true

  LeverageChanges.replicate(
    refreshcomp = (false || SOUP_2_NUTS),
    refreshtvolmcapcrsp = (false || SOUP_2_NUTS),
    refreshtretcrsp = (false || SOUP_2_NUTS),
    refreshdataseries=(false || SOUP_2_NUTS),
    refreshportfolios=(false || SOUP_2_NUTS),
    parallel = (true || SOUP_2_NUTS),
    trialrunonly = false, #WARNING MANUAL set to false on soup-2-nuts run

    #full set analysis tables
    table1 = (false || SOUP_2_NUTS),
    table2 = (false || SOUP_2_NUTS),
    table3 = (false || SOUP_2_NUTS),
    table4 = (false || SOUP_2_NUTS),
    table5 = (false || SOUP_2_NUTS),
    table6 = (false || SOUP_2_NUTS),

    #event analysis
    disttype=:primaryfull, #valid=[:primaryfull, :primaryshort, :oldfull, :oldshort]
    refreshseocrsp=(false  || SOUP_2_NUTS),
    refreshdistcrsp=(true  || SOUP_2_NUTS),
    refreshevent=(true  || SOUP_2_NUTS),
    refreshcumex=(true|| SOUP_2_NUTS),

    table7 = (true || SOUP_2_NUTS),
    table8 = (true || SOUP_2_NUTS),
    figure1 = (true  || SOUP_2_NUTS),
    figure2 = (true  || SOUP_2_NUTS),
    table9 = (true || SOUP_2_NUTS),
    figure3 = (true  || SOUP_2_NUTS),
    table10 = (true  || SOUP_2_NUTS),
    figure4 = (true  || SOUP_2_NUTS),

    tablea1 = (true || SOUP_2_NUTS),
    tablea2 = (true  || SOUP_2_NUTS),

    yearrange = 1961:2018,
    outsuffix = "all")

    println(LeverageChanges.YEAR_RANGE[])

end

#LeverageChanges.archivefile("lev-daily.jls.lz4")
#LeverageChanges.unarchivefile("lev-daily.jls.lz4")
#LeverageChanges.recompressfile("CCM")
#LeverageChanges.recompressfile("event\\crspdistributions")
#LeverageChanges.recompressfile("CRSP-D-0")
#LeverageChanges.recompressfile("CRSP-D-I")
#LeverageChanges.recompressfile("CRSP-D-II")
#LeverageChanges.recompressfile("COMP-A")
@time leveragechangesscript()

max_stars_repo_name: clintonTE/CCA

max_issues_repo_name: clintonTE/CCA

max_forks_repo_name: clintonTE/CCA


text:
function calc_R_0(cur_reactor::AbstractReactor)
  cur_constraint = cur_reactor.constraint

  cur_symbol = Symbol("R_0_from_$(cur_constraint)")

  cur_function = getfield(Fussy, cur_symbol)

  cur_R = cur_function(cur_reactor)

  cur_R
end


text:
"""
Line search method to apply once the direction is computed.
"""
abstract type LineSearchMethod end

"""
Computes step size: `2/(2 + t)` at iteration `t`.
"""
struct Agnostic <: LineSearchMethod end

"""
Computes step size via the
 [Backtracking Line Search](https://arxiv.org/pdf/1806.05123.pdf)
 method.
"""
struct Backtracking <: LineSearchMethod end

"""
Computes a step size via
[Golden Section Search](https://en.wikipedia.org/wiki/Golden-section_search).
"""
struct Goldenratio <: LineSearchMethod end

"""
Computes a step size for nonconvex functions: `1/sqrt(t + 1)`.
"""
struct Nonconvex <: LineSearchMethod end

"""
Computes the 'Short step' step size:
`dual_gap / (L * norm(x - v)^2)`,
where `L` is the Lipschitz constant of the gradient, `x` is the
current iterate, and `v` is the current Frank-Wolfe vertex.
"""
struct Shortstep <: LineSearchMethod end

"""
Constant step size given by `gamma0`
"""
struct FixedStep <: LineSearchMethod end

"""
Computes a 'Rational Short step' step size:
`sum((x - v) .* gradient ) // (L * sum((x - v) .^ 2))`,
where `L` is the Lipschitz constant of the gradient, `x` is the
current iterate, and `v` is the current Frank-Wolfe vertex.
"""
struct RationalShortstep <: LineSearchMethod end

"""
Slight modification of
Adaptive Step Size strategy from this
[paper](https://arxiv.org/pdf/1806.05123.pdf).
"""
struct Adaptive <: LineSearchMethod end

"""
Emphasis given to the algorithm for memory-saving or not.
The memory-saving mode may not be faster than the default
blas mode for small dimensions.
"""
@enum Emphasis blas = 1 memory = 2

@enum StepType begin
    initial = 1
    regular = 2
    lazy = 3
    lazylazy = 4
    dualstep = 5
    away = 6
    pairwise = 7
    drop = 8
    simplex_descent = 101
    gap_step = 102
    last = 1000
    pp = 1001
end

const st = (
    initial="I",
    regular="FW",
    lazy="L",
    lazylazy="LL",
    dualstep="LD",
    away="A",
    pairwise="P",
    drop="D",
    simplex_descent="SD",
    gap_step="GS",
    last="Last",
    pp="PP",
)


text:
export ft_setdatacharacteristics
export FT_BITS_8, FT_BITS_7, FT_STOP_BITS_1, FT_STOP_BITS_2
export FT_PARITY_NONE, FT_PARITY_ODD, FT_PARITY_EVEN, FT_PARITY_MARK
export FT_PARITY_SPACE

#Word Length (see ft_setdatacharacteristics)
const FT_BITS_8 = 8
const FT_BITS_7 = 7

#Stop Bits (see ft_setdatacharacteristics)
const FT_STOP_BITS_1 = 0
const FT_STOP_BITS_2 = 2

#Parity (see ft_setdatacharacteristics)
const FT_PARITY_NONE = 0
const FT_PARITY_ODD = 1
const FT_PARITY_EVEN = 2
const FT_PARITY_MARK = 3
const FT_PARITY_SPACE = 4

function ft_setdatacharacteristics(ft_handle::Culong,
                                   wordlength::Integer,
                                   stopbits::Integer,
                                   parity::Integer)
  @assert (wordlength==FT_BITS_8) || (wordlength==FT_BITS_7) "invalid word length"
  @assert (stopbits==FT_STOP_BITS_1) || (stopbits==FT_STOP_BITS_2) "invalid stop bits"
  @assert parity == FT_PARITY_NONE || parity == FT_PARITY_ODD ||
          parity == FT_PARITY_EVEN || parity == FT_PARITY_MARK ||
          parity == FT_PARITY_SPACE "invalid parity"
  ft_status = ccall((:FT_SetDataCharacteristics, d2xx),
                     Cuint,
                     (Culong, Cuchar, Cuchar, Cuchar),
                     ft_handle, wordlength, stopbits, parity)
  checkstatus(ft_status)
  return nothing
end


text:
module TestPipeline

using Random
using Test
using AMLPipelineBase
using AMLPipelineBase.Pipelines
using AMLPipelineBase.DecisionTreeLearners
using AMLPipelineBase.CrossValidators
using AMLPipelineBase.Utils

const data     = getiris()
const features = data[:,1:4]
const X        = data[:,1:5]
const Y        = data[:,5] |> Vector
X[!,5]         = X[!,5] .|> string

const ohe = OneHotEncoder()
const noop = Identity()
const rf  = RandomForest()
const ada = Adaboost()
const pt  = PrunedTree()
const numf  = NumFeatureSelector()
const catf  = CatFeatureSelector()

acc(X,Y) = score(:accuracy,X,Y)

function test_pipeline()
  # test initialization of types
  ohe = OneHotEncoder()
  linear1 = Pipeline(Dict(:name=>"lp",:machines => [ohe]))
  linear2 = Pipeline(Dict(:name=>"lp",:machines => [ohe]))
  combo1 = ComboPipeline(Dict(:machines=>[ohe,ohe]))
  combo2 = ComboPipeline(Dict(:machines=>[linear1,linear2]))
  linear1 = Pipeline([ohe])
  linear2 = Pipeline([ohe])
  combo1 = ComboPipeline([ohe,ohe])
  combo2 = ComboPipeline([linear1,linear2])
  # test fit/transform workflow
  res1 = fit_transform!(combo1,X)
  res2 = fit_transform!(combo2,X)
  @test (res1 .== res2) |> Matrix |> sum == 2100
  res1 = fit_transform(combo1,X)
  res2 = fit_transform(combo2,X)
  @test (res1 .== res2) |> Matrix |> sum == 2100
  res3=transform!(combo2,X)
  res4=fit_transform!(combo2,X)
  @test (res3 .== res4) |> Matrix |> sum == 2100
  pcombo1 = @pipeline ohe + ohe
  pres1 = fit_transform!(pcombo1,X)
  @test (pres1 .== res1) |> Matrix |> sum == 2100
  combo3 = ohe + ohe
  combo5 = linear1 + linear2
  res5 = fit_transform!(combo3,X)
  @test (res5 .== res2) |> Matrix |> sum == 2100
  res6 = fit_transform!(combo5,X)
  @test (res6 .== res2) |> Matrix |> sum == 2100
  res7 = fit_transform(combo5,X)
  @test (res6 .== res7) |> Matrix |> sum == 2100
end
@testset "Pipelines" begin
  Random.seed!(123)
  test_pipeline()
end

function test_sympipeline()
  pcombo5 = @pipeline :((ohe + noop) |> (ada * rf * pt))
  pcombo6 = (ohe + noop) |> (ada * rf * pt)
  @test crossvalidate(pcombo5,X,Y,acc,5,false).mean >= 0.90
  @test crossvalidate(pcombo6,X,Y,acc,5,false).mean >= 0.90
  expr = :((ohe + noop) |> (ada * rf * pt))
  processexpr!(expr.args)
  @test crossvalidate(eval(expr),X,Y,acc,5,false).mean >= 0.90
  pcombo6 = sympipeline(expr) |> eval
  @test crossvalidate(pcombo6,X,Y,acc,5,false).mean >= 0.90
  pcombo7 = (@pipelinez expr) |> eval
  @test crossvalidate(pcombo7,X,Y,acc,5,false).mean >= 0.90
end
@testset "Symbolic Pipeline: Global Scope" begin
  Random.seed!(123)
  test_sympipeline()
end

function test_advancedpipeline()
   expr = @pipeline ((((ada) + (rf) ) |> ohe) + numf) |> ada;
   @test crossvalidate(expr,X,Y,acc,5,false).mean >= 0.90
end
@testset "Advanced Pipeline: Learners as filters" begin
  Random.seed!(123)
  test_advancedpipeline()
end

function test_pipeline()
  # test symbolic pipeline expression
  pcombo2 = @pipeline ohe + noop
  pcombo3 = ohe + noop
  @test fit_transform!(pcombo2,features) |> Matrix |> size |> collect |> sum == 158
  @test fit_transform(pcombo2,features) |> Matrix |> size |> collect |> sum == 158
  @test fit_transform!(pcombo3,features) |> Matrix |> size |> collect |> sum == 158
  @test fit_transform(pcombo3,features) |> Matrix |> size |> collect |> sum == 158
  pcombo4 = @pipeline ohe + noop |> rf
  pcombo5 = (ohe + noop) |> rf
  @test crossvalidate(pcombo4,X,Y,acc,5,false).mean >= 0.90
  @test crossvalidate(pcombo5,X,Y,acc,5,false).mean >= 0.90
end
@testset "Symbolic Pipeline: Local Scope" begin
  Random.seed!(123)
  test_pipeline()
end

function test_split()
   dt1 = train_test_split(features,Y)
   perf1 = pipe_performance(rf,(x,y)->score(:accuracy,x,y), dt1.trX,dt1.trY,dt1.tstX,dt1.tstY)
   @test perf1 > 50.0
   dt2 = train_test_split(data[:,1:3],Vector(data[:,4]))
   perf2 = pipe_performance(rf,(x,y)->score(:rmse,x,y), dt2.trX,dt2.trY,dt2.tstX,dt2.tstY)
   @test perf2 < 0.50
end
@testset "Performance and train/test split" begin
  Random.seed!(123)
  test_split()
end

end


text:
@testset "time integration" begin
    D = [-1.0 2.0*π; -2.0*π -1.0]
    x0 = [1.0 0.0]
    test_data = (x0*exp(D))[:]
    e = ForwardEuler(1e-3)
    heuns = Heuns(1e-3)
    ssprk3 = StableRK3(1e-3)
    ssprk4 = StableRK4(1e-3)
    @test test_data≈integrate_time(x0,D,1.0,e) atol=e.step_size*10.0
    @test test_data≈integrate_time(x0,D,1.0,heuns) atol=heuns.step_size*10.0
    @test test_data≈integrate_time(x0,D,1.0,ssprk3) atol=ssprk3.step_size*10.0
    @test test_data≈integrate_time(x0,D,1.0,ssprk3) atol=ssprk4.step_size*10.0

    @testset "limiter -- MUSCL" begin
        # using DiscretisedFluidQueues
        one_D_T = [0.0][:,:]
        one_D_c = [1.0]
        one_D_model = BoundedFluidQueue(one_D_T,one_D_c,2.0)

        nodes = collect(0.0:2.0/50.0:2.0)

        mesh = DGMesh(nodes,7)#FRAPMesh(nodes,7)#

        dq = DiscretisedFluidQueue(one_D_model,mesh)
        B = build_full_generator(dq)
        # B[2,(end-n_bases_per_cell(mesh)+1):end] =
        #     B[n_bases_per_cell(mesh)+2,2:n_bases_per_cell(mesh)+1]
        # B[2,1] = 0.0
        t = 1.0
        f(x,i) = Float64(x<0.5)#sin(π*x)#
        d0 = SFMDistribution(f,dq)# interior_point_mass(-0.5,1,dq)#

        pdf0 = pdf(d0)
        pdf0_limit = pdf(limit(d0))
        x_vals = (-0.55:0.001:-0.45).+1
        @test any(pdf0.(x_vals,1).>1.01)
        @test !any(pdf0_limit.(x_vals,1).>1.01)
        @test any(pdf0.(x_vals,1).<-0.01)
        @test !any(pdf0_limit.(x_vals,1).<-0.01)

        dt = integrate_time(d0,B,t,StableRK4(0.01))
        dt_no_limit = integrate_time(d0,B,t,StableRK4(0.01); limiter=NoLimiter)

        pdft = pdf(dt)
        pdft_no_limit = pdf(limit(dt_no_limit))
        x_vals = (-0.05:0.001:0.55).+1
        @test !any(pdft.(x_vals,1).>1.01)
        @test any(pdft_no_limit.(x_vals,1).>1.01)
        @test !any(pdft.(x_vals,1).<-0.01)
        @test any(pdft_no_limit.(x_vals,1).<-0.01)
    end
end


text:
@testset "Small special tests" begin
@testset "Sum" begin
    com = CS.init()

    com_grid = Array{CS.Variable, 1}(undef, 7)
    com_grid[1] = add_var!(com, 1, 9)
    com_grid[2] = add_var!(com, 1, 9; fix=5)
    com_grid[3] = add_var!(com, 1, 9)
    com_grid[4] = add_var!(com, 1, 9)
    com_grid[5] = add_var!(com, 1, 9)

    com_grid[6] = add_var!(com, 1, 2)
    com_grid[7] = add_var!(com, 3, 5)

    CS.rm!(com, com_grid[4], 5)
    CS.remove_above!(com, com_grid[5], 2)

    add_constraint!(com, sum([com_grid[CartesianIndex(ind)] for ind in [1,2]]) == 11)
    add_constraint!(com, sum([com_grid[CartesianIndex(ind)] for ind in [3,4]]) == 11)
    add_constraint!(com, sum([com_grid[CartesianIndex(ind)] for ind in [2,5]]) ==  6)
    # testing coefficients from left and right
    add_constraint!(com, com_grid[6]*1+2*com_grid[7] ==  7)

    status = solve!(com; backtrack=false, keep_logs=true)
    # remove without pruning
    CS.remove_above!(com, com_grid[3], 5)

    # should also work for a fixed variable
    @test CS.compress_var_string(com_grid[1]) == "6"

    str_output = CS.get_str_repr(com_grid)
    @test str_output[1] == "6, 5, 2:5, [2, 3, 4, 6, 7, 8, 9], 1, 1, 3"

    com_grid_2D = [com_grid[1] com_grid[2]; com_grid[3] com_grid[4]]
    str_output = CS.get_str_repr(com_grid_2D)
    @test str_output[1] == "          6                    2:5          "
    @test str_output[2] == "          5           [2, 3, 4, 6, 7, 8, 9] "

    println(com_grid)

    @test status != :Infeasible
    @test CS.isfixed(com_grid[1])
    @test value(com_grid[1]) == 6
    @test CS.isfixed(com_grid[5])
    @test value(com_grid[5]) == 1
    @test !CS.has(com_grid[3], 6)

    @test CS.isfixed(com_grid[6])
    @test value(com_grid[6]) == 1

    @test CS.isfixed(com_grid[7])
    @test value(com_grid[7]) == 3
end

@testset "Reordering sum constraint" begin
    com = CS.init()

    x = add_var!(com, 0, 9)
    y = add_var!(com, 0, 9)
    z = add_var!(com, 0, 9)

    c1 = 2x+3x == 5
    @test length(c1.indices) == 1
    @test c1.indices[1] == 1
    @test c1.coeffs[1] == 5
    @test c1.rhs == 5

    add_constraint!(com, 2x+3x == 5)
    add_constraint!(com, 2x-3y+6+x == z)
    add_constraint!(com, x+2 == z)
    add_constraint!(com, z-2 == x)
    add_constraint!(com, 2x+x == z+3y-6)

    status = solve!(com)
    @test status == :Solved
    @test CS.isfixed(x) && value(x) == 1
    @test 2*value(x)-3*value(y)+6+value(x) == value(z)
    @test value(x)+2 == value(z)
end

@testset "Infeasible coeff sum" begin
    com = CS.init()

    com_grid = Array{CS.Variable, 1}(undef, 3)
    com_grid[1] = add_var!(com, 1, 9)
    com_grid[2] = add_var!(com, 1, 9)
    com_grid[3] = add_var!(com, 1, 9)

    CS.rm!(com, com_grid[2], 2)
    CS.rm!(com, com_grid[2], 3)
    CS.rm!(com, com_grid[2], 4)
    CS.rm!(com, com_grid[2], 6)
    CS.rm!(com, com_grid[2], 8)

    add_constraint!(com, com_grid[2]*2+5*com_grid[1]+0*com_grid[3] == 21)

    status = solve!(com; backtrack=true)
    @test status == :Infeasible
end

@testset "Negative coeff sum" begin
    com = CS.init()

    com_grid = Array{CS.Variable, 1}(undef, 3)
    com_grid[1] = add_var!(com, 1, 9)
    com_grid[2] = add_var!(com, 1, 9)
    com_grid[3] = add_var!(com, 1, 9)

    CS.remove_above!(com, com_grid[3], 3)

    add_constraint!(com, sum([7,5,-10].*com_grid) == -13)

    status = solve!(com; backtrack=true)
    @test status == :Solved
    @test CS.isfixed(com_grid[1]) && value(com_grid[1]) == 1
    @test CS.isfixed(com_grid[2]) && value(com_grid[2]) == 2
    @test CS.isfixed(com_grid[3]) && value(com_grid[3]) == 3

    com = CS.init()

    com_grid = Array{CS.Variable, 1}(undef, 3)
    com_grid[1] = add_var!(com, 1, 5)
    com_grid[2] = add_var!(com, 1, 2)
    com_grid[3] = add_var!(com, 1, 9)

    add_constraint!(com, sum([7,5,-10].*com_grid) == -13)

    status = solve!(com; backtrack=true)
    @test status == :Solved
    @test CS.isfixed(com_grid[1]) && value(com_grid[1]) == 1
    @test CS.isfixed(com_grid[2]) && value(com_grid[2]) == 2
    @test CS.isfixed(com_grid[3]) && value(com_grid[3]) == 3

    com = CS.init()

    v1 = add_var!(com, 1, 5)
    v2 = add_var!(com, 5, 10)

    add_constraint!(com, v1-v2 == 0)

    status = solve!(com; backtrack=false)
    @test status == :Solved
    @test CS.isfixed(v1) && value(v1) == 5
    @test CS.isfixed(v2) && value(v2) == 5
end

@testset "Equal constraint" begin
    # nothing to do
    com = CS.init()

    v1 = add_var!(com, 1, 2; fix=2)
    v2 = add_var!(com, 1, 2; fix=2)

    add_constraint!(com, v1 == v2)

    status = solve!(com)
    @test status == :Solved
    @test !com.info.backtracked
    @test CS.isfixed(v1) && value(v1) == 2
    @test CS.isfixed(v2) && value(v2) == 2

    # normal
    com = CS.init()

    v1 = add_var!(com, 1, 2)
    v2 = add_var!(com, 1, 2; fix=2)

    add_constraint!(com, v1 == v2)

    status = solve!(com)
    @test status == :Solved
    @test !com.info.backtracked
    @test CS.isfixed(v1) && value(v1) == 2
    @test CS.isfixed(v2) && value(v2) == 2

    # set but infeasible
    com = CS.init()

    v1 = add_var!(com, 1, 2)
    v2 = add_var!(com, 1, 2; fix=2)

    add_constraint!(com, v1 == v2)
    add_constraint!(com, CS.all_different([v1, v2]))

    status = solve!(com)
    @test status == :Infeasible
    @test !com.info.backtracked

    # set but infeasible reversed
    com = CS.init()

    v1 = add_var!(com, 1, 2)
    v2 = add_var!(com, 1, 2; fix=2)

    add_constraint!(com, v2 == v1)
    add_constraint!(com, CS.all_different([v1, v2]))

    status = solve!(com)
    @test status == :Infeasible
    @test !com.info.backtracked

    # reversed
    com = CS.init()

    v1 = add_var!(com, 1, 2)
    v2 = add_var!(com, 1, 2; fix=2)

    add_constraint!(com, v2 == v1)

    status = solve!(com)
    @test status == :Solved
    @test !com.info.backtracked
    @test CS.isfixed(v1) && value(v1) == 2
    @test CS.isfixed(v2) && value(v2) == 2

    # test with more
    com = CS.init()

    v1 = add_var!(com, 1, 2)
    v2 = add_var!(com, 1, 2; fix=2)
    v3 = add_var!(com, 1, 2)
    v4 = add_var!(com, 1, 2)

    add_constraint!(com, v1 == v2)
    add_constraint!(com, v1 == v4)
    add_constraint!(com, v1 == v3)

    status = solve!(com)
    @test status == :Solved
    @test !com.info.backtracked
    @test CS.isfixed(v1) && value(v1) == 2
    @test CS.isfixed(v2) && value(v2) == 2
    @test CS.isfixed(v3) && value(v3) == 2
    @test CS.isfixed(v4) && value(v4) == 2
end

@testset "NotSolved or infeasible" begin
    # NotSolved without backtracking
    com = CS.init()

    v1 = add_var!(com, 1, 2)
    v2 = add_var!(com, 1, 2)

    add_constraint!(com, v2 == v1)

    status = solve!(com; backtrack=false)
    @test status == :NotSolved
    @test !com.info.backtracked

    # Infeasible without backtracking
    com = CS.init()

    v1 = add_var!(com, 1, 2; fix=1)
    v2 = add_var!(com, 1, 2; fix=2)

    add_constraint!(com, v2 == v1)

    status = solve!(com)
    @test status == :Infeasible
    @test !com.info.backtracked

    # Infeasible without backtracking reverse
    com = CS.init()

    v1 = add_var!(com, 1, 2; fix=1)
    v2 = add_var!(com, 1, 2; fix=2)

    add_constraint!(com, v1 == v2)

    status = solve!(com)
    @test status == :Infeasible
    @test !com.info.backtracked
end

@testset "Test Equals()" begin
    # test using equal
    com = CS.init()

    v1 = add_var!(com, 1, 2)
    v2 = add_var!(com, 1, 2; fix=2)
    v3 = add_var!(com, 1, 2)
    v4 = add_var!(com, 1, 2)

    add_constraint!(com, CS.equal([v1,v2,v3,v4]))

    status = solve!(com)
    @test status == :Solved
    @test !com.info.backtracked
    @test CS.isfixed(v1) && value(v1) == 2
    @test CS.isfixed(v2) && value(v2) == 2
    @test CS.isfixed(v3) && value(v3) == 2
    @test CS.isfixed(v4) && value(v4) == 2
end

@testset "Test Equals() NotSolved/Infeasible" begin
    # test using equal
    com = CS.init()

    v1 = add_var!(com, 1, 2)
    v2 = add_var!(com, 1, 2; fix=2)
    v3 = add_var!(com, 1, 2)
    v4 = add_var!(com, 1, 3; fix=3)

    add_constraint!(com, CS.equal([v1,v2,v3,v4]))

    status = solve!(com)
    @test status == :Infeasible
    @test !com.info.backtracked
    @test CS.isfixed(v2) && value(v2) == 2
    @test CS.isfixed(v4) && value(v4) == 3

    # test using equal
    com = CS.init()

    v1 = add_var!(com, 1, 2)
    v2 = add_var!(com, 1, 2)
    v3 = add_var!(com, 1, 2)
    v4 = add_var!(com, 1, 2)

    add_constraint!(com, CS.equal([v1,v2,v3,v4]))

    status = solve!(com; backtrack=false)
    @test status == :NotSolved
    @test !com.info.backtracked
    @test !CS.isfixed(v2)
    @test !CS.isfixed(v4)
end

@testset "Graph coloring small" begin
    com = CS.init()

    # cover from numberphile video
    v1 = add_var!(com, 1, 4)
    v2 = add_var!(com, 1, 4)
    v3 = add_var!(com, 1, 4)
    v4 = add_var!(com, 1, 4)
    v5 = add_var!(com, 1, 4)
    v6 = add_var!(com, 1, 4)
    v7 = add_var!(com, 1, 4)
    v8 = add_var!(com, 1, 4)
    v9 = add_var!(com, 1, 4)

    add_constraint!(com, v1 != v2)
    add_constraint!(com, v1 != v3)
    add_constraint!(com, v1 != v4)
    add_constraint!(com, v1 != v5)
    add_constraint!(com, v2 != v3)
    add_constraint!(com, v2 != v5)
    add_constraint!(com, v4 != v3)
    add_constraint!(com, v4 != v5)
    add_constraint!(com, v2 != v6)
    add_constraint!(com, v2 != v7)
    add_constraint!(com, v3 != v7)
    add_constraint!(com, v3 != v8)
    add_constraint!(com, v4 != v8)
    add_constraint!(com, v4 != v9)
    add_constraint!(com, v5 != v9)
    add_constraint!(com, v5 != v6)
    add_constraint!(com, v5 != v6)
    add_constraint!(com, v6 != v9)
    add_constraint!(com, v6 != v9)
    add_constraint!(com, v6 != v7)
    add_constraint!(com, v8 != v9)
    add_constraint!(com, v8 != v7)

    status = solve!(com)
    @test status == :Solved
    @test com.info.backtracked

    # Infeasible with 3 colors
    com = CS.init()
    v1 = add_var!(com, 1, 3)
    v2 = add_var!(com, 1, 3)
    v3 = add_var!(com, 1, 3)
    v4 = add_var!(com, 1, 3)
    v5 = add_var!(com, 1, 3)
    v6 = add_var!(com, 1, 3)
    v7 = add_var!(com, 1, 3)
    v8 = add_var!(com, 1, 3)
    v9 = add_var!(com, 1, 3)

    add_constraint!(com, v1 != v2)
    add_constraint!(com, v1 != v3)
    add_constraint!(com, v1 != v4)
    add_constraint!(com, v1 != v5)
    add_constraint!(com, v2 != v3)
    add_constraint!(com, v2 != v5)
    add_constraint!(com, v4 != v3)
    add_constraint!(com, v4 != v5)
    add_constraint!(com, v2 != v6)
    add_constraint!(com, v2 != v7)
    add_constraint!(com, v3 != v7)
    add_constraint!(com, v3 != v8)
    add_constraint!(com, v4 != v8)
    add_constraint!(com, v4 != v9)
    add_constraint!(com, v5 != v9)
    add_constraint!(com, v5 != v6)
    add_constraint!(com, v5 != v6)
    add_constraint!(com, v6 != v9)
    add_constraint!(com, v6 != v9)
    add_constraint!(com, v6 != v7)
    add_constraint!(com, v8 != v9)
    add_constraint!(com, v8 != v7)

    status = solve!(com)
    @test status == :Infeasible
    @test com.info.backtracked

    # very simple
    com = CS.init()

    v1 = add_var!(com, 1, 2; fix=1)
    v2 = add_var!(com, 1, 2)

    add_constraint!(com, v1 != v2)

    status = solve!(com; backtrack=false)
    @test status == :Solved
    @test !com.info.backtracked
    @test CS.isfixed(v1) && value(v1) == 1
    @test CS.isfixed(v2) && value(v2) == 2

    # infeasible at beginning
    com = CS.init()

    v1 = add_var!(com, 1, 2; fix=1)
    v2 = add_var!(com, 1, 2; fix=1)

    add_constraint!(com, v1 != v2)

    status = solve!(com; backtrack=false)
    @test status == :Infeasible
    @test !com.info.backtracked

    # Error not implemented
    com = CS.init()

    v1 = add_var!(com, 1, 2; fix=1)
    v2 = add_var!(com, 1, 2)
    v3 = add_var!(com, 1, 2)

    @test_throws ErrorException add_constraint!(com, !CS.equal([v1,v2,v3]))
    @test_throws ErrorException add_constraint!(com, !CS.all_different([v1,v2,v3]))
end

end


text:
using LambdaTAME: ΛTAME_param_search_profiled, ΛTAME_param_search,  ΛTAME

@testset "ΛTAME" begin

    β = 0.0
    α = 1.0
    tol = 1e-6
    max_iter = 15

    @testset "Type Stability" begin

        @suppress_out begin
            # --  Third Order Symmetric Tensors  -- #
            @inferred ΛTAME_param_search_profiled(A_TOST,B_TOST)
            @inferred ΛTAME_param_search(A_TOST,B_TOST)
            @inferred ΛTAME(A_TOST,  B_TOST,β,max_iter,tol,α)

            @inferred ΛTAME_param_search_profiled(A_UST,B_UST)
            @inferred ΛTAME_param_search(A_UST,B_UST)
            @inferred ΛTAME(A_UST,  B_UST,β,max_iter,tol,α)

            # --  Unweighted Symmetric Tensors Multi Motif -- #
            @inferred ΛTAME_param_search(A_UST_MM,B_UST_MM)
            @inferred ΛTAME(A_UST_MM,  B_UST_MM,β,max_iter,tol,α)
        end
    end

    @testset "Types output comparisons" begin

        # shared

        #=
        # LowRankTAME + TAME only
        d = 5
        U = rand(A_TOST.n,d)
        V = rand(B_TOST.n,d)
        X = U*V'
        =#
        @suppress_out begin
            #no shifts
            U_TOST,  V_TOST   = ΛTAME(A_TOST,  B_TOST,β,max_iter,tol,α)
            U_UST,   V_UST    = ΛTAME(A_UST,   B_UST ,β,max_iter,tol,α)
            U_UST_MM,V_UST_MM = ΛTAME([A_UST],[B_UST],β,max_iter,tol,α) #multimotif routines

            @test norm(U_TOST - U_UST)/norm(U_UST) < round_off_bound(U_TOST)
            @test norm(U_UST_MM - U_UST)/norm(U_UST_MM) < round_off_bound(U_UST_MM)

            @test norm(V_TOST - V_UST)/norm(V_UST) < round_off_bound(V_TOST)
            @test norm(V_UST_MM - V_UST)/norm(V_UST_MM)< round_off_bound(V_UST_MM)

            #shifts
            β = 1.0
            α = .5
            U_TOST,  V_TOST   = ΛTAME(A_TOST,  B_TOST,β,max_iter,tol,α)
            U_UST,   V_UST    = ΛTAME(A_UST,   B_UST ,β,max_iter,tol,α)
            U_UST_MM,V_UST_MM = ΛTAME([A_UST],[B_UST],β,max_iter,tol,α) #multimotif routines

            @test norm(U_TOST - U_UST)/norm(U_UST) < round_off_bound(U_TOST)
            @test norm(U_UST_MM - U_UST)/norm(U_UST_MM) < round_off_bound(U_UST_MM)

            @test norm(V_TOST - V_UST)/norm(V_UST) < round_off_bound(V_TOST)
            @test norm(V_UST_MM - V_UST)/norm(V_UST_MM)< round_off_bound(V_UST_MM)

        end
    end

end

max_stars_repo_name: chuckcol/LambdaTAME

max_issues_repo_name: chuckcol/LambdaTAME

max_forks_repo_name: chuckcol/LambdaTAME