text:
using LibGit2

"""
    download(; repo = default_repo(), path = default_path())

Download a workshop from github.

Keyword arguments `repo` and `path` can be used to specify the github repository
and the local installation directory.
"""
function download(
    ;
    repo = default_repo(),
    path = default_path(),
    auto_overwrite = false,
)
    if !occursin("github.com", repo)
        repo = joinpath(GITHUB_BASEURL, repo)
    end
    workshop = basename(repo)
    target = joinpath(path, workshop)

    # overwrite?
    if isdir(target)
        if auto_overwrite
            answer = true
        else
            @info "Target directory $(target) already exists. Should I overwrite?"
            answer = yes_no_dialog()
        end

        if answer == true
            rm(target, force = true, recursive = true)
        else
            @info "Aborting."
            return false
        end
    end

    @info "Downloading \"$workshop\" to $target"
    LibGit2.clone(repo, target)
    return true
end

"""
    install(; repo = default_repo(), path = default_path(), global_IJulia = true, auto_overwrite = true)

More programmatic workshop installation interface.

By default, the workshop will be downloaded to the desktop (on windows)
or the home directory (on linux/macOS). The installation path
can be adjusted per keyword argument `path = desired/install/path`.
"""
function install(
    ;
    repo = default_repo(),
    path = default_path(),
    global_IJulia = true,
    auto_overwrite = true,
)
    success = download(
        repo = repo,
        path = path,
        auto_overwrite = auto_overwrite,
    )
    !success && (return false)
    install_dependencies(joinpath(path, basename(repo)))
    workshop_dir = joinpath(path, basename(repo))
    @info "Workshop successfully installed to $(workshop_dir)."
    println()

    if global_IJulia
        install_IJulia()
    end
    return true
end

"""
    install_dependencies(workshop_path)

Instantiates and precompiles the workshops Julia environment (Manifest.toml).
Depending on the workshop and your machine this may take a while (even >15 minutes).
"""
function install_dependencies(workshop_path)
    @info "Installing and precompiling all dependencies (this may take a while)"
    with_pkg_env(workshop_path, change_dir = true) do
        pkg"instantiate"
        pkg"precompile"
        println("Done.")
        println()
    end
end

"""
    run_wizard()

Start the workshop installation wizard. It will interactively guide you through
the selection, download, and installation of a Julia workshop.

If not already present, the wizard will also install IJulia.
"""
function run_wizard()
    @info "Welcome to the workshop installation wizard."

    # default settings?
    @info "Should I download and install the $(default_workshop()) workshop with default settings?"
    answer = yes_no_dialog()
    if answer == true
        success = install(
            repo = default_repo(),
            path = default_path(),
            global_IJulia = true,
            auto_overwrite = false,
        )
    else
        success = _install_interactive()
    end

    !success && (return nothing)
    @info "That's it. Start the notebook server with \"using IJulia; notebook()\" and take it away."
    println()
    println("- The Workshop Wizard")
    return nothing
end

function _install_interactive()
    # workshop selection
    @info "Please select a workshop:"
    workshop = workshop_selector()
    if isnothing(workshop)
        @info "No workshop selected. Aborting."
        return false
    end
    @info "Selected \"$workshop\" for installation."

    if !julia_version_is_compatible(workshop)
        ver = WORKSHOP_JULIA_VERSION[workshop]
        @warn "Your Julia version ($(VERSION)) does not match the recommended version " *
              "for this workshop ($(ver))!\n" *
              "THE WORKSHOP MAY NOT INSTALL/RUN PROPERLY!\n" *
              "To fix this, upgrade/downgrade to Julia version $(ver)."
    end

    # installation directory
    path = default_path()
    @info "The default installation directory is: $path"
    @info "Is that okay?"
    answer = yes_no_dialog()
    if answer == false
        @info "Please specify a directory:"
        p = ""
        while !isdir(p)
            p = readline()
            if !isdir(p)
                @info "Invalid directory. Please try again."
            end
        end
        path = p
    end

    @info "Starting the installation"
    success = install(
        repo = joinpath(GITHUB_BASEURL, workshop),
        path = path,
        global_IJulia = false,
        auto_overwrite = false,
    )
    !success && (return false)

    @info "Should I install IJulia globally for you?"
    answer = yes_no_dialog()
    if answer == true
        install_IJulia()
    end
    return true
end


text:
const CACHE = 256
buf = zeros(Int,CACHE)
buf[1] = 1
#happy(n) returns 1 if happy, 0 if not
function happy(n)
	if n < CACHE
		buf[n] > 0 && return 2-buf[n]
		buf[n] = 2
	end
	sum = 0
	nn = n
	while nn != 0
		x = nn%10
		sum += x*x
		nn = int8(nn/10)
	end
	x = happy(sum)
	n < CACHE && (buf[n] = 2-x)
	return x
end
function main()
	i = 1; counter = 1000000
	while counter > 0
		if happy(i) == 1
			counter -= 1
		end
		i += 1
	end
	return i-1
end

ext: jl

max_stars_repo_name: ethansaxenian/RosettaDecode

max_issues_repo_name: ethansaxenian/RosettaDecode

max_forks_repo_name: ethansaxenian/RosettaDecode


text:
module ArbRadixFloatingPoints

using LinearAlgebra
using StaticArrays

export ArbRadixFloat

import Base: +, -, *, /, \, ^, abs, conj, convert, float, imag, precision, real
import LinearAlgebra: norm

"""
Floating point with arbitrary radix

# parameters
radix : "base" of floating point number
precision: number of digits
Tdigit : type of digit. If not specified in convert(), defaults to Int

significand and exponent are stored as fields
"""
struct ArbRadixFloat{radix,precision,Tdigit} <: Number
    significand :: SVector{precision,Tdigit}
    exponent :: Int
end

include("convert.jl")
include("show.jl")

# Basic floating point interface

precision(num::ArbRadixFloat{radix,precision,Tdigit}) where {radix,precision,Tdigit} = precision

conj(num::ArbRadixFloat{radix,precision,Tdigit}) where {radix,precision,Tdigit} =
    ArbRadixFloat{conj(radix),precision,Tdigit}(num.significand, num.exponent)

-(num::ArbRadixFloat{radix,precision,Tdigit}) where {radix,precision,Tdigit} =
    ArbRadixFloat{conj(radix),precision,Tdigit}(-num.significand, num.exponent)

function float(num::ArbRadixFloat{radix,precision,Tdigit}) where radix where precision where Tdigit
    anstype = promote_type(typeof(float(radix)), Tdigit)
    convert(anstype, ArbRadixFloat{float(radix),precision,Tdigit}(num.significand, num.exponent))
end

# These operations are just placeholders in lieu of actually doing the computations natively

for func in (:abs, :imag, :real, :norm)
    @eval begin
        function ($func)(num::ArbRadixFloat{radix,precision,Tdigit}) where {radix,precision,Tdigit}
            $(func)(convert(promote_type(typeof(radix), Tdigit), num))
        end
    end
end

for op in (:+, :-, :*, :/, :\, :^)
    @eval begin
        function ($op)(x::ArbRadixFloat{radix,precision,Tdigit}, y::ArbRadixFloat{radix,precision,Tdigit}) where {radix,precision,Tdigit}
            Twork = promote_type(typeof(radix), Tdigit)
            convert(ArbRadixFloat{radix,precision,Tdigit}, ($op)(convert(Twork, x), convert(Twork, y)))
        end
    end
end

end # module


text:
module SLEEF

export xsin,  xcos,  xtan,  xasin,  xacos,  xatan,  xlog,
	xexp,  xsinh,  xcosh,  xtanh,  xasinh,  xacosh,  xatanh,
	xcbrt,  xexp2,  xexp10,  xexpm1,  xlog10,  xlog1p,
	xsin_u1, xcos_u1,  xtan_u1,  xasin_u1,
	xacos_u1,  xatan_u1,  xlog_u1,  xcbrt_u1,
	xsincos, xsincos_u1, xatan2, xatan2_u1, xpow,
	xldexp, xilogb

if isfile(joinpath(dirname(@__FILE__),"..","deps","deps.jl"))
    include("../deps/deps.jl")
else
    error("SLEEF not properly installed. Please run Pkg.build(\"SLEEF\")")
end

include("sleef_inner.jl")

##############################

end # SLEEF Module


text:
load("dag/pshow_dag.jl")

c=CallNode(SymNode(:+,:call), LiteralNode(1), LiteralNode(2))
d=CallNode(SymNode(:+,:call), c,c)

pprintln(d)

max_stars_repo_name: toivoh/julia-kernels

max_issues_repo_name: toivoh/julia-kernels

max_forks_repo_name: toivoh/julia-kernels


text:
##### Arithmetics on contiguous ranks #####

for m=0:4, n=0:4
    global addrank
    @eval addrank(::Type{ContRank{$m}}, ::Type{ContRank{$n}}) = ContRank{$(m+n)}
end

addrank(::Type{ContRank{M}}, ::Type{ContRank{N}}) where {M,N} = ContRank{M+N}
addrank(::Type{ContRank{0}}, ::Type{ContRank{N}}) where {N} = ContRank{N}
addrank(::Type{ContRank{N}}, ::Type{ContRank{0}}) where {N} = ContRank{N}

for m=0:4, n=0:4
    global minrank
    @eval minrank(::Type{ContRank{$m}}, ::Type{ContRank{$n}}) = ContRank{$(min(m,n))}
end

minrank(::Type{ContRank{M}}, ::Type{ContRank{N}}) where {M,N} = ContRank{min(M,N)}
minrank(::Type{ContRank{0}}, ::Type{ContRank{N}}) where {N} = ContRank{0}
minrank(::Type{ContRank{N}}, ::Type{ContRank{0}}) where {N} = ContRank{0}

for m=0:4, n=0:4
    global restrict_crank
    @eval restrict_crank(::Type{ContRank{$m}}, ::NTuple{$n,Int}) = ContRank{$(min(m,n))}
end

restrict_crank(::Type{ContRank{M}}, ::NTuple{N,Int}) where {M,N} = ContRank{min(M,N)}
restrict_crank(::Type{ContRank{0}}, ::NTuple{N,Int}) where {N} = ContRank{0}
restrict_crank(::Type{ContRank{N}}, ::NTuple{0,Int}) where {N} = ContRank{0}

### contiguous rank computation based on indices

# 0D
contrank() = ContRank{0}

# 1D
contrank(i::Union{Colon,UnitRange}) = ContRank{1}
contrank(i::Subs) = ContRank{0}

# 2D
contrank(i1::Colon, i2::Union{Colon,UnitRange}) = ContRank{2}
contrank(i1::Colon, i2::Subs) = ContRank{1}
contrank(i1::UnitRange, i2::Subs) = ContRank{1}
contrank(i1::Subs, i2::Subs) = ContRank{0}

# 3D
contrank(i1::Colon, i2::Colon, i3::Colon) = ContRank{3}
contrank(i1::Colon, i2::Colon, i3::UnitRange) = ContRank{3}
contrank(i1::Colon, i2::Colon, i3::Subs) = ContRank{2}
contrank(i1::Colon, i2::UnitRange, i3::Subs) = ContRank{2}
contrank(i1::Colon, i2::Subs, i3::Subs) = ContRank{1}
contrank(i1::UnitRange, i2::Subs, i3::Subs) = ContRank{1}
contrank(i1::Subs, i2::Subs, i3::Subs) = ContRank{0}

# ND
contrank(i1::Colon, i2::Colon, i3::Colon, I::Subs...) = addrank(ContRank{3}, contrank(I...))
contrank(i1::Colon, i2::Colon, i3::UnitRange, I::Subs...) = ContRank{3}
contrank(i1::Colon, i2::Colon, i3::Subs, I::Subs...) = ContRank{2}
contrank(i1::Colon, i2::UnitRange, i3::Subs, I::Subs...) = ContRank{2}
contrank(i1::Colon, i2::Subs, i3::Subs, I::Subs...) = ContRank{1}
contrank(i1::UnitRange, i2::Subs, i3::Subs, I::Subs...) = ContRank{1}
contrank(i1::Subs, i2::Subs, i3::Subs, I::Subs...) = ContRank{0}

# contiguous rank with array & arrayviews

acontrank(a::Array, i1::Subs) = contrank(i1)
acontrank(a::Array, i1::Subs, i2::Subs) = contrank(i1, i2)
acontrank(a::Array, i1::Subs, i2::Subs, i3::Subs) = contrank(i1, i2, i3)
acontrank(a::Array, i1::Subs, i2::Subs, i3::Subs, i4::Subs) = contrank(i1, i2, i3, i4)
acontrank(a::Array, i1::Subs, i2::Subs, i3::Subs, i4::Subs, i5::Subs, I::Subs...) =
    contrank(i1, i2, i3, i4, i5, I...)

acontrank(a::StridedArrayView{T,N,N}, i1::Subs) where {T,N} = contrank(i1)
acontrank(a::StridedArrayView{T,N,N}, i1::Subs, i2::Subs) where {T,N} = contrank(i1, i2)
acontrank(a::StridedArrayView{T,N,N}, i1::Subs, i2::Subs, i3::Subs) where {T,N} = contrank(i1, i2, i3)
acontrank(a::StridedArrayView{T,N,N}, i1::Subs, i2::Subs, i3::Subs, i4::Subs) where {T,N} =
    contrank(i1, i2, i3, i4)
acontrank(a::StridedArrayView{T,N,N}, i1::Subs, i2::Subs, i3::Subs, i4::Subs, i5::Subs, I::Subs...) where {T,N} =
    contrank(i1, i2, i3, i4, i5, I...)

acontrank(a::StridedArrayView{T,N,M}, i1::Subs) where {T,N,M} = minrank(contrank(i1), ContRank{M})
acontrank(a::StridedArrayView{T,N,M}, i1::Subs, i2::Subs) where {T,N,M} =
    minrank(contrank(i1, i2), ContRank{M})
acontrank(a::StridedArrayView{T,N,M}, i1::Subs, i2::Subs, i3::Subs) where {T,N,M} =
    minrank(contrank(i1, i2, i3), ContRank{M})
acontrank(a::StridedArrayView{T,N,M}, i1::Subs, i2::Subs, i3::Subs, i4::Subs) where {T,N,M} =
    minrank(contrank(i1, i2, i3, i4), ContRank{M})
acontrank(a::StridedArrayView{T,N,M}, i1::Subs, i2::Subs, i3::Subs, i4::Subs, i5::Subs, I::Subs...) where {T,N,M} =
    minrank(contrank(i1, i2, i3, i4, i5, I...), ContRank{M})


text:
# Autogenerated wrapper script for LibModbus_jll for aarch64-linux-gnu
export libmodbus

JLLWrappers.@generate_wrapper_header("LibModbus")
JLLWrappers.@declare_library_product(libmodbus, "libmodbus.so.5")
function __init__()
    JLLWrappers.@generate_init_header()
    JLLWrappers.@init_library_product(
        libmodbus,
        "lib/libmodbus.so",
        RTLD_LAZY | RTLD_DEEPBIND,
    )

    JLLWrappers.@generate_init_footer()
end  # __init__()


text:
using Rubin
using Tests
using Elliptic
using HypergeometricFunctions
using Polylogarithms
using SpecialFunctions

@test integrate(x^7*fresnels(b*x), x) == :(1//8*x^8*fresnels(b*x)+-105//8*pi^-4*b^-8*fresnels(b*x)+-35//8*pi^-3*b^-5*x^3*cos((1/2)*pi*b^2*x^2)+-7//8*pi^-2*b^-3*x^5*sin((1/2)*pi*b^2*x^2)+1//8*pi^-1*b^-1*x^7*cos((1/2)*pi*b^2*x^2)+105//8*x*pi^-4*b^-7*sin((1/2)*pi*b^2*x^2))
@test integrate(x^6*fresnels(b*x), x) == :(1//7*x^7*fresnels(b*x)+48//7*pi^-4*b^-7*sin((1/2)*pi*b^2*x^2)+-24//7*pi^-3*b^-5*x^2*cos((1/2)*pi*b^2*x^2)+-6//7*pi^-2*b^-3*x^4*sin((1/2)*pi*b^2*x^2)+1//7*pi^-1*b^-1*x^6*cos((1/2)*pi*b^2*x^2))
@test integrate(x^5*fresnels(b*x), x) == :(1//6*x^6*fresnels(b*x)+5//2*pi^-3*b^-6*fresnelc(b*x)+-5//2*x*pi^-3*b^-5*cos((1/2)*pi*b^2*x^2)+-5//6*pi^-2*b^-3*x^3*sin((1/2)*pi*b^2*x^2)+1//6*pi^-1*b^-1*x^5*cos((1/2)*pi*b^2*x^2))
@test integrate(x^4*fresnels(b*x), x) == :(1//5*x^5*fresnels(b*x)+-8//5*pi^-3*b^-5*cos((1/2)*pi*b^2*x^2)+-4//5*pi^-2*b^-3*x^2*sin((1/2)*pi*b^2*x^2)+1//5*pi^-1*b^-1*x^4*cos((1/2)*pi*b^2*x^2))
@test integrate(x^3*fresnels(b*x), x) == :(1//4*x^4*fresnels(b*x)+3//4*pi^-2*b^-4*fresnels(b*x)+-3//4*x*pi^-2*b^-3*sin((1/2)*pi*b^2*x^2)+1//4*pi^-1*b^-1*x^3*cos((1/2)*pi*b^2*x^2))
@test integrate(x^2*fresnels(b*x), x) == :(1//3*x^3*fresnels(b*x)+-2//3*pi^-2*b^-3*sin((1/2)*pi*b^2*x^2)+1//3*pi^-1*b^-1*x^2*cos((1/2)*pi*b^2*x^2))
@test integrate(x*fresnels(b*x), x) == :((1/2)*x^2*fresnels(b*x)+-1//2*pi^-1*b^-2*fresnelc(b*x)+(1/2)*x*pi^-1*b^-1*cos((1/2)*pi*b^2*x^2))
@test integrate(fresnels(b*x), x) == :(x*fresnels(b*x)+pi^-1*b^-1*cos((1/2)*pi*b^2*x^2))
@test integrate(x^-2*fresnels(b*x), x) == :((1/2)*b*sinintegral((1/2)*pi*b^2*x^2)+-1*x^-1*fresnels(b*x))
@test integrate(x^-3*fresnels(b*x), x) == :(-1//2*x^-2*fresnels(b*x)+(1/2)*pi*b^2*fresnelc(b*x)+-1//2*b*x^-1*sin((1/2)*pi*b^2*x^2))
@test integrate(x^-4*fresnels(b*x), x) == :(-1//3*x^-3*fresnels(b*x)+-1//6*b*x^-2*sin((1/2)*pi*b^2*x^2)+1//12*pi*b^3*cosintegral((1/2)*pi*b^2*x^2))
@test integrate(x^-5*fresnels(b*x), x) == :(-1//4*x^-4*fresnels(b*x)+-1//12*b*x^-3*sin((1/2)*pi*b^2*x^2)+-1//12*pi^2*b^4*fresnels(b*x)+-1//12*pi*b^3*x^-1*cos((1/2)*pi*b^2*x^2))
@test integrate(x^-6*fresnels(b*x), x) == :(-1//5*x^-5*fresnels(b*x)+-1//20*b*x^-4*sin((1/2)*pi*b^2*x^2)+-1//80*pi^2*b^5*sinintegral((1/2)*pi*b^2*x^2)+-1//40*pi*b^3*x^-2*cos((1/2)*pi*b^2*x^2))
@test integrate(x^-7*fresnels(b*x), x) == :(-1//6*x^-6*fresnels(b*x)+-1//30*b*x^-5*sin((1/2)*pi*b^2*x^2)+-1//90*pi^3*b^6*fresnelc(b*x)+-1//90*pi*b^3*x^-3*cos((1/2)*pi*b^2*x^2)+1//90*pi^2*b^5*x^-1*sin((1/2)*pi*b^2*x^2))
@test integrate(x^-8*fresnels(b*x), x) == :(-1//7*x^-7*fresnels(b*x)+-1//42*b*x^-6*sin((1/2)*pi*b^2*x^2)+-1//672*pi^3*b^7*cosintegral((1/2)*pi*b^2*x^2)+-1//168*pi*b^3*x^-4*cos((1/2)*pi*b^2*x^2)+1//336*pi^2*b^5*x^-2*sin((1/2)*pi*b^2*x^2))
@test integrate(x^-9*fresnels(b*x), x) == :(-1//8*x^-8*fresnels(b*x)+-1//56*b*x^-7*sin((1/2)*pi*b^2*x^2)+1//840*pi^4*b^8*fresnels(b*x)+-1//280*pi*b^3*x^-5*cos((1/2)*pi*b^2*x^2)+1//840*pi^2*b^5*x^-3*sin((1/2)*pi*b^2*x^2)+1//840*pi^3*b^7*x^-1*cos((1/2)*pi*b^2*x^2))
@test integrate(x^-10*fresnels(b*x), x) == :(-1//9*x^-9*fresnels(b*x)+-1//72*b*x^-8*sin((1/2)*pi*b^2*x^2)+1//6912*pi^4*b^9*sinintegral((1/2)*pi*b^2*x^2)+-1//432*pi*b^3*x^-6*cos((1/2)*pi*b^2*x^2)+1//1728*pi^2*b^5*x^-4*sin((1/2)*pi*b^2*x^2)+1//3456*pi^3*b^7*x^-2*cos((1/2)*pi*b^2*x^2))
@test integrate((c+d*x)^3*fresnels(a+b*x), x) == :(1//4*d^-1*(c+d*x)^4*fresnels(a+b*x)+pi^-1*b^-4*(b*c+-1*a*d)^3*cos((1/2)*pi*(a+b*x)^2)+-1//4*b^-4*d^-1*(b*c+-1*a*d)^4*fresnels(a+b*x)+3//4*pi^-2*b^-4*d^3*fresnels(a+b*x)+-2*pi^-2*b^-4*d^2*(b*c+-1*a*d)*sin((1/2)*pi*(a+b*x)^2)+-3//2*d*pi^-1*b^-4*(b*c+-1*a*d)^2*fresnelc(a+b*x)+-3//4*pi^-2*b^-4*d^3*(a+b*x)*sin((1/2)*pi*(a+b*x)^2)+1//4*pi^-1*b^-4*d^3*(a+b*x)^3*cos((1/2)*pi*(a+b*x)^2)+pi^-1*b^-4*d^2*(a+b*x)^2*(b*c+-1*a*d)*cos((1/2)*pi*(a+b*x)^2)+3//2*d*pi^-1*b^-4*(b*c+-1*a*d)^2*(a+b*x)*cos((1/2)*pi*(a+b*x)^2))
@test integrate((c+d*x)^2*fresnels(a+b*x), x) == :(1//3*d^-1*(c+d*x)^3*fresnels(a+b*x)+pi^-1*b^-3*(b*c+-1*a*d)^2*cos((1/2)*pi*(a+b*x)^2)+-2//3*pi^-2*b^-3*d^2*sin((1/2)*pi*(a+b*x)^2)+-1//3*b^-3*d^-1*(b*c+-1*a*d)^3*fresnels(a+b*x)+-1*d*pi^-1*b^-3*(b*c+-1*a*d)*fresnelc(a+b*x)+1//3*pi^-1*b^-3*d^2*(a+b*x)^2*cos((1/2)*pi*(a+b*x)^2)+d*pi^-1*b^-3*(a+b*x)*(b*c+-1*a*d)*cos((1/2)*pi*(a+b*x)^2))
@test integrate((c+d*x)*fresnels(a+b*x), x) == :((1/2)*d^-1*(c+d*x)^2*fresnels(a+b*x)+pi^-1*b^-2*(b*c+-1*a*d)*cos((1/2)*pi*(a+b*x)^2)+-1//2*d*pi^-1*b^-2*fresnelc(a+b*x)+-1//2*b^-2*d^-1*(b*c+-1*a*d)^2*fresnels(a+b*x)+(1/2)*d*pi^-1*b^-2*(a+b*x)*cos((1/2)*pi*(a+b*x)^2))
@test integrate(fresnels(a+b*x), x) == :(pi^-1*b^-1*cos((1/2)*pi*(a+b*x)^2)+b^-1*(a+b*x)*fresnels(a+b*x))
@test integrate(x^3*fresnels(a+b*x), x) == :(1//4*x^4*fresnels(a+b*x)+-1//4*a^4*b^-4*fresnels(a+b*x)+3//4*pi^-2*b^-4*fresnels(a+b*x)+-1*pi^-1*a^3*b^-4*cos((1/2)*pi*(a+b*x)^2)+2*a*pi^-2*b^-4*sin((1/2)*pi*(a+b*x)^2)+-3//2*pi^-1*a^2*b^-4*fresnelc(a+b*x)+-1//4*pi^-2*b^-4*(3a+3*b*x)*sin((1/2)*pi*(a+b*x)^2)+1//4*pi^-1*b^-4*(a+b*x)^3*cos((1/2)*pi*(a+b*x)^2)+-1*a*pi^-1*b^-4*(a+b*x)^2*cos((1/2)*pi*(a+b*x)^2)+3//2*pi^-1*a^2*b^-4*(a+b*x)*cos((1/2)*pi*(a+b*x)^2))
@test integrate(x^2*fresnels(a+b*x), x) == :(1//3*x^3*fresnels(a+b*x)+-2//3*pi^-2*b^-3*sin((1/2)*pi*(a+b*x)^2)+1//3*a^3*b^-3*fresnels(a+b*x)+a*pi^-1*b^-3*fresnelc(a+b*x)+pi^-1*a^2*b^-3*cos((1/2)*pi*(a+b*x)^2)+1//3*pi^-1*b^-3*(a+b*x)^2*cos((1/2)*pi*(a+b*x)^2)+-1*a*pi^-1*b^-3*(a+b*x)*cos((1/2)*pi*(a+b*x)^2))
@test integrate(x*fresnels(a+b*x), x) == :((1/2)*x^2*fresnels(a+b*x)+-1//2*pi^-1*b^-2*fresnelc(a+b*x)+-1//2*a^2*b^-2*fresnels(a+b*x)+(1/2)*pi^-1*b^-2*(a+b*x)*cos((1/2)*pi*(a+b*x)^2)+-1*a*pi^-1*b^-2*cos((1/2)*pi*(a+b*x)^2))
@test integrate(fresnels(a+b*x), x) == :(pi^-1*b^-1*cos((1/2)*pi*(a+b*x)^2)+b^-1*(a+b*x)*fresnels(a+b*x))
@test integrate(x^7*fresnels(b*x)^2, x) == :(1//8*x^8*fresnels(b*x)^2+10*pi^-5*b^-8*sin(pi*b^2*x^2)+-105//8*pi^-4*b^-8*fresnels(b*x)^2+-105//16*pi^-4*b^-6*x^2+7//48*pi^-2*b^-2*x^6+-55//16*pi^-4*b^-6*x^2*cos(pi*b^2*x^2)+-5//8*pi^-3*b^-4*x^4*sin(pi*b^2*x^2)+1//16*pi^-2*b^-2*x^6*cos(pi*b^2*x^2)+-35//4*pi^-3*b^-5*x^3*cos((1/2)*pi*b^2*x^2)*fresnels(b*x)+-7//4*pi^-2*b^-3*x^5*fresnels(b*x)*sin((1/2)*pi*b^2*x^2)+1//4*pi^-1*b^-1*x^7*cos((1/2)*pi*b^2*x^2)*fresnels(b*x)+105//4*x*pi^-4*b^-7*fresnels(b*x)*sin((1/2)*pi*b^2*x^2))
@test integrate(x^6*fresnels(b*x)^2, x) == :(1//7*x^7*fresnels(b*x)^2+-48//7*x*pi^-4*b^-6+6//35*pi^-2*b^-2*x^5+-21//8*x*pi^-4*b^-6*cos(pi*b^2*x^2)+-17//28*pi^-3*b^-4*x^3*sin(pi*b^2*x^2)+1//14*pi^-2*b^-2*x^5*cos(pi*b^2*x^2)+96//7*pi^-4*b^-7*fresnels(b*x)*sin((1/2)*pi*b^2*x^2)+531//112*2^(1/2)*pi^-4*b^-7*fresnelc(b*x*2^(1/2))+-48//7*pi^-3*b^-5*x^2*cos((1/2)*pi*b^2*x^2)*fresnels(b*x)+-12//7*pi^-2*b^-3*x^4*fresnels(b*x)*sin((1/2)*pi*b^2*x^2)+2//7*pi^-1*b^-1*x^6*cos((1/2)*pi*b^2*x^2)*fresnels(b*x))
@test integrate(x^4*fresnels(b*x)^2, x) == :(1//5*x^5*fresnels(b*x)^2+4//15*pi^-2*b^-2*x^3+-16//5*pi^-3*b^-5*cos((1/2)*pi*b^2*x^2)*fresnels(b*x)+-11//20*x*pi^-3*b^-4*sin(pi*b^2*x^2)+1//10*pi^-2*b^-2*x^3*cos(pi*b^2*x^2)+43//40*2^(1/2)*pi^-3*b^-5*fresnels(b*x*2^(1/2))+-8//5*pi^-2*b^-3*x^2*fresnels(b*x)*sin((1/2)*pi*b^2*x^2)+2//5*pi^-1*b^-1*x^4*cos((1/2)*pi*b^2*x^2)*fresnels(b*x))
@test integrate(x^3*fresnels(b*x)^2, x) == :(1//4*x^4*fresnels(b*x)^2+-1//2*pi^-3*b^-4*sin(pi*b^2*x^2)+3//4*pi^-2*b^-4*fresnels(b*x)^2+3//8*pi^-2*b^-2*x^2+1//8*pi^-2*b^-2*x^2*cos(pi*b^2*x^2)+(1/2)*pi^-1*b^-1*x^3*cos((1/2)*pi*b^2*x^2)*fresnels(b*x)+-3//2*x*pi^-2*b^-3*fresnels(b*x)*sin((1/2)*pi*b^2*x^2))
@test integrate(x^2*fresnels(b*x)^2, x) == :(1//3*x^3*fresnels(b*x)^2+2//3*x*pi^-2*b^-2+-5//12*2^(1/2)*pi^-2*b^-3*fresnelc(b*x*2^(1/2))+-4//3*pi^-2*b^-3*fresnels(b*x)*sin((1/2)*pi*b^2*x^2)+1//6*x*pi^-2*b^-2*cos(pi*b^2*x^2)+2//3*pi^-1*b^-1*x^2*cos((1/2)*pi*b^2*x^2)*fresnels(b*x))
@test integrate(fresnels(b*x)^2, x) == :(x*fresnels(b*x)^2+2*pi^-1*b^-1*cos((1/2)*pi*b^2*x^2)*fresnels(b*x)+-1//2*2^(1/2)*pi^-1*b^-1*fresnels(b*x*2^(1/2)))
@test integrate(x^-5*fresnels(b*x)^2, x) == :(-1//4*x^-4*fresnels(b*x)^2+-1//24*b^2*x^-2+-1//12*pi^2*b^4*fresnels(b*x)^2+1//12*pi*b^4*sinintegral(pi*b^2*x^2)+1//24*b^2*x^-2*cos(pi*b^2*x^2)+-1//6*b*x^-3*fresnels(b*x)*sin((1/2)*pi*b^2*x^2)+-1//6*pi*b^3*x^-1*cos((1/2)*pi*b^2*x^2)*fresnels(b*x))
@test integrate(x^-9*fresnels(b*x)^2, x) == :(-1//8*x^-8*fresnels(b*x)^2+-1//336*b^2*x^-6+-1//280*pi^3*b^8*sinintegral(pi*b^2*x^2)+1//336*b^2*x^-6*cos(pi*b^2*x^2)+1//840*pi^4*b^8*fresnels(b*x)^2+1//1680*pi^2*b^6*x^-2+-1//28*b*x^-7*fresnels(b*x)*sin((1/2)*pi*b^2*x^2)+-1//336*pi^2*b^6*x^-2*cos(pi*b^2*x^2)+-1//420*pi*b^4*x^-4*sin(pi*b^2*x^2)+-1//140*pi*b^3*x^-5*cos((1/2)*pi*b^2*x^2)*fresnels(b*x)+1//420*pi^2*b^5*x^-3*fresnels(b*x)*sin((1/2)*pi*b^2*x^2)+1//420*pi^3*b^7*x^-1*cos((1/2)*pi*b^2*x^2)*fresnels(b*x))
@test integrate(fresnels(a+b*x)^2, x) == :(b^-1*fresnels(a+b*x)^2*(a+b*x)+2*pi^-1*b^-1*cos((1/2)*pi*(a+b*x)^2)*fresnels(a+b*x)+-1//2*2^(1/2)*pi^-1*b^-1*fresnels(2^(1/2)*(a+b*x)))
@test integrate(x^2*fresnels(d*(a+b*log(c*x^n))), x) == :(1//3*x^3*fresnels(d*(a+b*log(c*x^n)))+x^3*(c*x^n)^(-3*n^-1)*(1//12+-1//12*im)*erf(pi^-1//2*b^-1*d^-1*(1/2+(1/2)*im)*(3*n^-1+pi*im*a*b*d^2+pi*im*b^2*d^2*log(c*x^n)))*exp(-3*a*b^-1*n^-1+9//2*im*pi^-1*b^-2*d^-2*n^-2)+x^3*(c*x^n)^(-3*n^-1)*(1//12+-1//12*im)*erfi(pi^-1//2*b^-1*d^-1*(1/2+(1/2)*im)*(3*n^-1+-1*pi*im*a*b*d^2+-1*pi*im*b^2*d^2*log(c*x^n)))*exp(-3*a*b^-1*n^-1+-9//2*im*pi^-1*b^-2*d^-2*n^-2))
@test integrate(x*fresnels(d*(a+b*log(c*x^n))), x) == :((1/2)*x^2*fresnels(d*(a+b*log(c*x^n)))+x^2*(c*x^n)^(-2*n^-1)*(1//8+-1//8*im)*erf(pi^-1//2*b^-1*d^-1*(1/2+(1/2)*im)*(2*n^-1+pi*im*a*b*d^2+pi*im*b^2*d^2*log(c*x^n)))*exp(pi^-1*b^-2*d^-2*n^-2*(2im+-2*pi*a*b*n*d^2))+x^2*(c*x^n)^(-2*n^-1)*(1//8+-1//8*im)*erfi(pi^-1//2*b^-1*d^-1*(1/2+(1/2)*im)*(2*n^-1+-1*pi*im*a*b*d^2+-1*pi*im*b^2*d^2*log(c*x^n)))*exp(-1*pi^-1*b^-2*d^-2*n^-2*(2im+2*pi*a*b*n*d^2)))
@test integrate(fresnels(d*(a+b*log(c*x^n))), x) == :(x*fresnels(d*(a+b*log(c*x^n)))+x*(c*x^n)^(-1*n^-1)*(1//4+-1//4*im)*erf(pi^-1//2*b^-1*d^-1*(1/2+(1/2)*im)*(n^-1+pi*im*a*b*d^2+pi*im*b^2*d^2*log(c*x^n)))*exp(-1//2*b^-2*n^-2*(-1*im*pi^-1*d^-2+2*a*b*n))+x*(c*x^n)^(-1*n^-1)*(1//4+-1//4*im)*erfi(pi^-1//2*b^-1*d^-1*(1/2+(1/2)*im)*(n^-1+-1*pi*im*a*b*d^2+-1*pi*im*b^2*d^2*log(c*x^n)))*exp(-1//2*b^-2*n^-2*(im*pi^-1*d^-2+2*a*b*n)))
@test integrate(x^-1*fresnels(d*(a+b*log(c*x^n))), x) == :(b^-1*n^-1*(a+b*log(c*x^n))*fresnels(d*(a+b*log(c*x^n)))+pi^-1*b^-1*d^-1*n^-1*cos((1/2)*pi*d^2*(a+b*log(c*x^n))^2))
@test integrate(x^-2*fresnels(d*(a+b*log(c*x^n))), x) == :(-1*x^-1*fresnels(d*(a+b*log(c*x^n)))+x^-1*(c*x^n)^(n^-1)*(1//4+-1//4*im)*erf(pi^-1//2*b^-1*d^-1*(1/2+(1/2)*im)*(n^-1+-1*pi*im*a*b*d^2+-1*pi*im*b^2*d^2*log(c*x^n)))*exp((1/2)*b^-2*n^-2*(im*pi^-1*d^-2+2*a*b*n))+x^-1*(c*x^n)^(n^-1)*(1//4+-1//4*im)*erfi(pi^-1//2*b^-1*d^-1*(1/2+(1/2)*im)*(n^-1+pi*im*a*b*d^2+pi*im*b^2*d^2*log(c*x^n)))*exp((1/2)*b^-2*n^-2*(-1*im*pi^-1*d^-2+2*a*b*n)))
@test integrate(x^-3*fresnels(d*(a+b*log(c*x^n))), x) == :(-1//2*x^-2*fresnels(d*(a+b*log(c*x^n)))+x^-2*(c*x^n)^(2*n^-1)*(1//8+-1//8*im)*erf(pi^-1//2*b^-1*d^-1*(1/2+(1/2)*im)*(2*n^-1+-1*pi*im*a*b*d^2+-1*pi*im*b^2*d^2*log(c*x^n)))*exp(pi^-1*b^-2*d^-2*n^-2*(2im+2*pi*a*b*n*d^2))+x^-2*(c*x^n)^(2*n^-1)*(1//8+-1//8*im)*erfi(pi^-1//2*b^-1*d^-1*(1/2+(1/2)*im)*(2*n^-1+pi*im*a*b*d^2+pi*im*b^2*d^2*log(c*x^n)))*exp(-1*pi^-1*b^-2*d^-2*n^-2*(2im+-2*pi*a*b*n*d^2)))
@test integrate((e*x)^m*fresnels(d*(a+b*log(c*x^n))), x) == :(e^-1*(e*x)^(1+m)*(1+m)^-1*fresnels(d*(a+b*log(c*x^n)))+x*(c*x^n)^(-1*n^-1*(1+m))*(e*x)^m*(1+m)^-1*(1//4+-1//4*im)*erf(pi^-1//2*b^-1*d^-1*n^-1*(1/2+(1/2)*im)*(1+m+pi*im*a*b*n*d^2+pi*im*n*b^2*d^2*log(c*x^n)))*exp((1/2)*im*pi^-1*b^-2*d^-2*n^-2*(1+m)*(1+m+2*pi*im*a*b*n*d^2))+x*(c*x^n)^(-1*n^-1*(1+m))*(e*x)^m*(1+m)^-1*(1//4+-1//4*im)*erfi(pi^-1//2*b^-1*d^-1*n^-1*(1/2+(1/2)*im)*(1+m+-1*pi*im*a*b*n*d^2+-1*pi*im*n*b^2*d^2*log(c*x^n)))*exp(-1//2*im*pi^-1*b^-2*d^-2*n^-2*(1+m)*(1+m+-2*pi*im*a*b*n*d^2)))
@test integrate(fresnels(b*x)^2*sin((1/2)*pi*b^2*x^2), x) == :(1//3*b^-1*fresnels(b*x)^3)
@test integrate(fresnels(b*x)*sin((1/2)*pi*b^2*x^2), x) == :((1/2)*b^-1*fresnels(b*x)^2)
@test integrate(fresnels(b*x)^-1*sin((1/2)*pi*b^2*x^2), x) == :(b^-1*log(fresnels(b*x)))
@test integrate(fresnels(b*x)^-2*sin((1/2)*pi*b^2*x^2), x) == :(-1*b^-1*fresnels(b*x)^-1)
@test integrate(fresnels(b*x)^-3*sin((1/2)*pi*b^2*x^2), x) == :(-1//2*b^-1*fresnels(b*x)^-2)
@test integrate(fresnels(b*x)^n*sin((1/2)*pi*b^2*x^2), x) == :(b^-1*(1+n)^-1*fresnels(b*x)^(1+n))
@test integrate(x^8*fresnels(b*x)*sin((1/2)*pi*b^2*x^2), x) == :(-40*pi^-5*b^-9*sin(pi*b^2*x^2)+-7//12*pi^-2*b^-3*x^6+105//2*pi^-4*b^-9*fresnels(b*x)^2+105//4*pi^-4*b^-7*x^2+-1//4*pi^-2*b^-3*x^6*cos(pi*b^2*x^2)+5//2*pi^-3*b^-5*x^4*sin(pi*b^2*x^2)+55//4*pi^-4*b^-7*x^2*cos(pi*b^2*x^2)+-1*pi^-1*b^-2*x^7*cos((1/2)*pi*b^2*x^2)*fresnels(b*x)+-105*x*pi^-4*b^-8*fresnels(b*x)*sin((1/2)*pi*b^2*x^2)+7*pi^-2*b^-4*x^5*fresnels(b*x)*sin((1/2)*pi*b^2*x^2)+35*pi^-3*b^-6*x^3*cos((1/2)*pi*b^2*x^2)*fresnels(b*x))
@test integrate(x^7*fresnels(b*x)*sin((1/2)*pi*b^2*x^2), x) == :(24*x*pi^-4*b^-7+-3//5*pi^-2*b^-3*x^5+-48*pi^-4*b^-8*fresnels(b*x)*sin((1/2)*pi*b^2*x^2)+-531//32*2^(1/2)*pi^-4*b^-8*fresnelc(b*x*2^(1/2))+-1//4*pi^-2*b^-3*x^5*cos(pi*b^2*x^2)+17//8*pi^-3*b^-5*x^3*sin(pi*b^2*x^2)+147//16*x*pi^-4*b^-7*cos(pi*b^2*x^2)+-1*pi^-1*b^-2*x^6*cos((1/2)*pi*b^2*x^2)*fresnels(b*x)+6*pi^-2*b^-4*x^4*fresnels(b*x)*sin((1/2)*pi*b^2*x^2)+24*pi^-3*b^-6*x^2*cos((1/2)*pi*b^2*x^2)*fresnels(b*x))
@test integrate(x^5*fresnels(b*x)*sin((1/2)*pi*b^2*x^2), x) == :(-2//3*pi^-2*b^-3*x^3+8*pi^-3*b^-6*cos((1/2)*pi*b^2*x^2)*fresnels(b*x)+-43//16*2^(1/2)*pi^-3*b^-6*fresnels(b*x*2^(1/2))+-1//4*pi^-2*b^-3*x^3*cos(pi*b^2*x^2)+11//8*x*pi^-3*b^-5*sin(pi*b^2*x^2)+-1*pi^-1*b^-2*x^4*cos((1/2)*pi*b^2*x^2)*fresnels(b*x)+4*pi^-2*b^-4*x^2*fresnels(b*x)*sin((1/2)*pi*b^2*x^2))
@test integrate(x^4*fresnels(b*x)*sin((1/2)*pi*b^2*x^2), x) == :(pi^-3*b^-5*sin(pi*b^2*x^2)+-3//2*pi^-2*b^-5*fresnels(b*x)^2+-3//4*pi^-2*b^-3*x^2+-1//4*pi^-2*b^-3*x^2*cos(pi*b^2*x^2)+-1*pi^-1*b^-2*x^3*cos((1/2)*pi*b^2*x^2)*fresnels(b*x)+3*x*pi^-2*b^-4*fresnels(b*x)*sin((1/2)*pi*b^2*x^2))
@test integrate(x^3*fresnels(b*x)*sin((1/2)*pi*b^2*x^2), x) == :(-1*x*pi^-2*b^-3+2*pi^-2*b^-4*fresnels(b*x)*sin((1/2)*pi*b^2*x^2)+-1//4*x*pi^-2*b^-3*cos(pi*b^2*x^2)+5//8*2^(1/2)*pi^-2*b^-4*fresnelc(b*x*2^(1/2))+-1*pi^-1*b^-2*x^2*cos((1/2)*pi*b^2*x^2)*fresnels(b*x))
@test integrate(x*fresnels(b*x)*sin((1/2)*pi*b^2*x^2), x) == :(-1*pi^-1*b^-2*cos((1/2)*pi*b^2*x^2)*fresnels(b*x)+1//4*2^(1/2)*pi^-1*b^-2*fresnels(b*x*2^(1/2)))
@test integrate(fresnels(b*x)*sin((1/2)*pi*b^2*x^2), x) == :((1/2)*b^-1*fresnels(b*x)^2)
@test integrate(x^-4*fresnels(b*x)*sin((1/2)*pi*b^2*x^2), x) == :(-1//12*b*x^-2+-1//3*x^-3*fresnels(b*x)*sin((1/2)*pi*b^2*x^2)+-1//6*pi^2*b^3*fresnels(b*x)^2+1//6*pi*b^3*sinintegral(pi*b^2*x^2)+1//12*b*x^-2*cos(pi*b^2*x^2)+-1//3*pi*b^2*x^-1*cos((1/2)*pi*b^2*x^2)*fresnels(b*x))
@test integrate(x^-8*fresnels(b*x)*sin((1/2)*pi*b^2*x^2), x) == :(-1//84*b*x^-6+-1//7*x^-7*fresnels(b*x)*sin((1/2)*pi*b^2*x^2)+-1//70*pi^3*b^7*sinintegral(pi*b^2*x^2)+1//84*b*x^-6*cos(pi*b^2*x^2)+1//210*pi^4*b^7*fresnels(b*x)^2+1//420*pi^2*b^5*x^-2+-1//84*pi^2*b^5*x^-2*cos(pi*b^2*x^2)+-1//105*pi*b^3*x^-4*sin(pi*b^2*x^2)+-1//35*pi*b^2*x^-5*cos((1/2)*pi*b^2*x^2)*fresnels(b*x)+1//105*pi^2*b^4*x^-3*fresnels(b*x)*sin((1/2)*pi*b^2*x^2)+1//105*pi^3*b^6*x^-1*cos((1/2)*pi*b^2*x^2)*fresnels(b*x))
@test integrate(x^7*cos((1/2)*pi*b^2*x^2)*fresnels(b*x), x) == :(4*pi^-3*b^-5*x^3+-1//14*pi^-1*b^-1*x^7+-48*pi^-4*b^-8*cos((1/2)*pi*b^2*x^2)*fresnels(b*x)+-147//16*x*pi^-4*b^-7*sin(pi*b^2*x^2)+1//4*pi^-2*b^-3*x^5*sin(pi*b^2*x^2)+17//8*pi^-3*b^-5*x^3*cos(pi*b^2*x^2)+531//32*2^(1/2)*pi^-4*b^-8*fresnels(b*x*2^(1/2))+pi^-1*b^-2*x^6*fresnels(b*x)*sin((1/2)*pi*b^2*x^2)+-24*pi^-3*b^-6*x^2*fresnels(b*x)*sin((1/2)*pi*b^2*x^2)+6*pi^-2*b^-4*x^4*cos((1/2)*pi*b^2*x^2)*fresnels(b*x))
@test integrate(x^6*cos((1/2)*pi*b^2*x^2)*fresnels(b*x), x) == :(-11//2*pi^-4*b^-7*sin(pi*b^2*x^2)+-1//12*pi^-1*b^-1*x^6+15//2*pi^-3*b^-7*fresnels(b*x)^2+15//4*pi^-3*b^-5*x^2+1//4*pi^-2*b^-3*x^4*sin(pi*b^2*x^2)+7//4*pi^-3*b^-5*x^2*cos(pi*b^2*x^2)+pi^-1*b^-2*x^5*fresnels(b*x)*sin((1/2)*pi*b^2*x^2)+-15*x*pi^-3*b^-6*fresnels(b*x)*sin((1/2)*pi*b^2*x^2)+5*pi^-2*b^-4*x^3*cos((1/2)*pi*b^2*x^2)*fresnels(b*x))
@test integrate(x^5*cos((1/2)*pi*b^2*x^2)*fresnels(b*x), x) == :(4*x*pi^-3*b^-5+-1//10*pi^-1*b^-1*x^5+-8*pi^-3*b^-6*fresnels(b*x)*sin((1/2)*pi*b^2*x^2)+-43//16*2^(1/2)*pi^-3*b^-6*fresnelc(b*x*2^(1/2))+1//4*pi^-2*b^-3*x^3*sin(pi*b^2*x^2)+11//8*x*pi^-3*b^-5*cos(pi*b^2*x^2)+pi^-1*b^-2*x^4*fresnels(b*x)*sin((1/2)*pi*b^2*x^2)+4*pi^-2*b^-4*x^2*cos((1/2)*pi*b^2*x^2)*fresnels(b*x))
@test integrate(x^3*cos((1/2)*pi*b^2*x^2)*fresnels(b*x), x) == :(-1//6*pi^-1*b^-1*x^3+2*pi^-2*b^-4*cos((1/2)*pi*b^2*x^2)*fresnels(b*x)+-5//8*2^(1/2)*pi^-2*b^-4*fresnels(b*x*2^(1/2))+1//4*x*pi^-2*b^-3*sin(pi*b^2*x^2)+pi^-1*b^-2*x^2*fresnels(b*x)*sin((1/2)*pi*b^2*x^2))
@test integrate(x^2*cos((1/2)*pi*b^2*x^2)*fresnels(b*x), x) == :(-1//2*pi^-1*b^-3*fresnels(b*x)^2+-1//4*pi^-1*b^-1*x^2+1//4*pi^-2*b^-3*sin(pi*b^2*x^2)+x*pi^-1*b^-2*fresnels(b*x)*sin((1/2)*pi*b^2*x^2))
@test integrate(x*cos((1/2)*pi*b^2*x^2)*fresnels(b*x), x) == :(-1//2*x*pi^-1*b^-1+pi^-1*b^-2*fresnels(b*x)*sin((1/2)*pi*b^2*x^2)+1//4*2^(1/2)*pi^-1*b^-2*fresnelc(b*x*2^(1/2)))
@test integrate(x^-2*cos((1/2)*pi*b^2*x^2)*fresnels(b*x), x) == :(1//4*b*sinintegral(pi*b^2*x^2)+-1*x^-1*cos((1/2)*pi*b^2*x^2)*fresnels(b*x)+-1//2*pi*b*fresnels(b*x)^2)
@test integrate(x^-6*cos((1/2)*pi*b^2*x^2)*fresnels(b*x), x) == :(-7//120*pi^2*b^5*sinintegral(pi*b^2*x^2)+-1//5*x^-5*cos((1/2)*pi*b^2*x^2)*fresnels(b*x)+-1//40*b*x^-4*sin(pi*b^2*x^2)+1//30*pi^3*b^5*fresnels(b*x)^2+1//60*pi*b^3*x^-2+-1//24*pi*b^3*x^-2*cos(pi*b^2*x^2)+1//15*pi*b^2*x^-3*fresnels(b*x)*sin((1/2)*pi*b^2*x^2)+1//15*pi^2*b^4*x^-1*cos((1/2)*pi*b^2*x^2)*fresnels(b*x))
@test integrate(x^-10*cos((1/2)*pi*b^2*x^2)*fresnels(b*x), x) == :(-1//9*x^-9*cos((1/2)*pi*b^2*x^2)*fresnels(b*x)+-1//144*b*x^-8*sin(pi*b^2*x^2)+-1//1890*pi^5*b^9*fresnels(b*x)^2+-1//3780*pi^3*b^7*x^-2+1//756*pi*b^3*x^-6+83//30240*pi^4*b^9*sinintegral(pi*b^2*x^2)+-11//3024*pi*b^3*x^-6*cos(pi*b^2*x^2)+5//2016*pi^3*b^7*x^-2*cos(pi*b^2*x^2)+67//30240*pi^2*b^5*x^-4*sin(pi*b^2*x^2)+-1//945*pi^3*b^6*x^-3*fresnels(b*x)*sin((1/2)*pi*b^2*x^2)+-1//945*pi^4*b^8*x^-1*cos((1/2)*pi*b^2*x^2)*fresnels(b*x)+1//63*pi*b^2*x^-7*fresnels(b*x)*sin((1/2)*pi*b^2*x^2)+1//315*pi^2*b^4*x^-5*cos((1/2)*pi*b^2*x^2)*fresnels(b*x))
@test integrate(x^7*fresnelc(b*x), x) == :(1//8*x^8*fresnelc(b*x)+-105//8*pi^-4*b^-8*fresnelc(b*x)+-7//8*pi^-2*b^-3*x^5*cos((1/2)*pi*b^2*x^2)+-1//8*pi^-1*b^-1*x^7*sin((1/2)*pi*b^2*x^2)+35//8*pi^-3*b^-5*x^3*sin((1/2)*pi*b^2*x^2)+105//8*x*pi^-4*b^-7*cos((1/2)*pi*b^2*x^2))
@test integrate(x^6*fresnelc(b*x), x) == :(1//7*x^7*fresnelc(b*x)+48//7*pi^-4*b^-7*cos((1/2)*pi*b^2*x^2)+-6//7*pi^-2*b^-3*x^4*cos((1/2)*pi*b^2*x^2)+-1//7*pi^-1*b^-1*x^6*sin((1/2)*pi*b^2*x^2)+24//7*pi^-3*b^-5*x^2*sin((1/2)*pi*b^2*x^2))
@test integrate(x^5*fresnelc(b*x), x) == :(1//6*x^6*fresnelc(b*x)+-5//2*pi^-3*b^-6*fresnels(b*x)+-5//6*pi^-2*b^-3*x^3*cos((1/2)*pi*b^2*x^2)+-1//6*pi^-1*b^-1*x^5*sin((1/2)*pi*b^2*x^2)+5//2*x*pi^-3*b^-5*sin((1/2)*pi*b^2*x^2))
@test integrate(x^4*fresnelc(b*x), x) == :(1//5*x^5*fresnelc(b*x)+8//5*pi^-3*b^-5*sin((1/2)*pi*b^2*x^2)+-4//5*pi^-2*b^-3*x^2*cos((1/2)*pi*b^2*x^2)+-1//5*pi^-1*b^-1*x^4*sin((1/2)*pi*b^2*x^2))
@test integrate(x^3*fresnelc(b*x), x) == :(1//4*x^4*fresnelc(b*x)+3//4*pi^-2*b^-4*fresnelc(b*x)+-3//4*x*pi^-2*b^-3*cos((1/2)*pi*b^2*x^2)+-1//4*pi^-1*b^-1*x^3*sin((1/2)*pi*b^2*x^2))
@test integrate(x^2*fresnelc(b*x), x) == :(1//3*x^3*fresnelc(b*x)+-2//3*pi^-2*b^-3*cos((1/2)*pi*b^2*x^2)+-1//3*pi^-1*b^-1*x^2*sin((1/2)*pi*b^2*x^2))
@test integrate(x*fresnelc(b*x), x) == :((1/2)*x^2*fresnelc(b*x)+(1/2)*pi^-1*b^-2*fresnels(b*x)+-1//2*x*pi^-1*b^-1*sin((1/2)*pi*b^2*x^2))
@test integrate(fresnelc(b*x), x) == :(x*fresnelc(b*x)+-1*pi^-1*b^-1*sin((1/2)*pi*b^2*x^2))
@test integrate(x^-2*fresnelc(b*x), x) == :((1/2)*b*cosintegral((1/2)*pi*b^2*x^2)+-1*x^-1*fresnelc(b*x))
@test integrate(x^-3*fresnelc(b*x), x) == :(-1//2*x^-2*fresnelc(b*x)+-1//2*pi*b^2*fresnels(b*x)+-1//2*b*x^-1*cos((1/2)*pi*b^2*x^2))
@test integrate(x^-4*fresnelc(b*x), x) == :(-1//3*x^-3*fresnelc(b*x)+-1//6*b*x^-2*cos((1/2)*pi*b^2*x^2)+-1//12*pi*b^3*sinintegral((1/2)*pi*b^2*x^2))
@test integrate(x^-5*fresnelc(b*x), x) == :(-1//4*x^-4*fresnelc(b*x)+-1//12*b*x^-3*cos((1/2)*pi*b^2*x^2)+-1//12*pi^2*b^4*fresnelc(b*x)+1//12*pi*b^3*x^-1*sin((1/2)*pi*b^2*x^2))
@test integrate(x^-6*fresnelc(b*x), x) == :(-1//5*x^-5*fresnelc(b*x)+-1//20*b*x^-4*cos((1/2)*pi*b^2*x^2)+-1//80*pi^2*b^5*cosintegral((1/2)*pi*b^2*x^2)+1//40*pi*b^3*x^-2*sin((1/2)*pi*b^2*x^2))
@test integrate(x^-7*fresnelc(b*x), x) == :(-1//6*x^-6*fresnelc(b*x)+-1//30*b*x^-5*cos((1/2)*pi*b^2*x^2)+1//90*pi^3*b^6*fresnels(b*x)+1//90*pi*b^3*x^-3*sin((1/2)*pi*b^2*x^2)+1//90*pi^2*b^5*x^-1*cos((1/2)*pi*b^2*x^2))
@test integrate(x^-8*fresnelc(b*x), x) == :(-1//7*x^-7*fresnelc(b*x)+-1//42*b*x^-6*cos((1/2)*pi*b^2*x^2)+1//672*pi^3*b^7*sinintegral((1/2)*pi*b^2*x^2)+1//168*pi*b^3*x^-4*sin((1/2)*pi*b^2*x^2)+1//336*pi^2*b^5*x^-2*cos((1/2)*pi*b^2*x^2))
@test integrate(x^-9*fresnelc(b*x), x) == :(-1//8*x^-8*fresnelc(b*x)+-1//56*b*x^-7*cos((1/2)*pi*b^2*x^2)+1//840*pi^4*b^8*fresnelc(b*x)+-1//840*pi^3*b^7*x^-1*sin((1/2)*pi*b^2*x^2)+1//280*pi*b^3*x^-5*sin((1/2)*pi*b^2*x^2)+1//840*pi^2*b^5*x^-3*cos((1/2)*pi*b^2*x^2))
@test integrate(x^-10*fresnelc(b*x), x) == :(-1//9*x^-9*fresnelc(b*x)+-1//72*b*x^-8*cos((1/2)*pi*b^2*x^2)+1//6912*pi^4*b^9*cosintegral((1/2)*pi*b^2*x^2)+-1//3456*pi^3*b^7*x^-2*sin((1/2)*pi*b^2*x^2)+1//432*pi*b^3*x^-6*sin((1/2)*pi*b^2*x^2)+1//1728*pi^2*b^5*x^-4*cos((1/2)*pi*b^2*x^2))
@test integrate((c+d*x)^3*fresnelc(a+b*x), x) == :(1//4*d^-1*(c+d*x)^4*fresnelc(a+b*x)+-1*pi^-1*b^-4*(b*c+-1*a*d)^3*sin((1/2)*pi*(a+b*x)^2)+-1//4*b^-4*d^-1*(b*c+-1*a*d)^4*fresnelc(a+b*x)+3//4*pi^-2*b^-4*d^3*fresnelc(a+b*x)+-2*pi^-2*b^-4*d^2*(b*c+-1*a*d)*cos((1/2)*pi*(a+b*x)^2)+-3//4*pi^-2*b^-4*d^3*(a+b*x)*cos((1/2)*pi*(a+b*x)^2)+-1//4*pi^-1*b^-4*d^3*(a+b*x)^3*sin((1/2)*pi*(a+b*x)^2)+3//2*d*pi^-1*b^-4*(b*c+-1*a*d)^2*fresnels(a+b*x)+-1*pi^-1*b^-4*d^2*(a+b*x)^2*(b*c+-1*a*d)*sin((1/2)*pi*(a+b*x)^2)+-3//2*d*pi^-1*b^-4*(b*c+-1*a*d)^2*(a+b*x)*sin((1/2)*pi*(a+b*x)^2))
@test integrate((c+d*x)^2*fresnelc(a+b*x), x) == :(1//3*d^-1*(c+d*x)^3*fresnelc(a+b*x)+-1*pi^-1*b^-3*(b*c+-1*a*d)^2*sin((1/2)*pi*(a+b*x)^2)+-2//3*pi^-2*b^-3*d^2*cos((1/2)*pi*(a+b*x)^2)+-1//3*b^-3*d^-1*(b*c+-1*a*d)^3*fresnelc(a+b*x)+d*pi^-1*b^-3*(b*c+-1*a*d)*fresnels(a+b*x)+-1//3*pi^-1*b^-3*d^2*(a+b*x)^2*sin((1/2)*pi*(a+b*x)^2)+-1*d*pi^-1*b^-3*(a+b*x)*(b*c+-1*a*d)*sin((1/2)*pi*(a+b*x)^2))
@test integrate((c+d*x)*fresnelc(a+b*x), x) == :((1/2)*d^-1*(c+d*x)^2*fresnelc(a+b*x)+(1/2)*d*pi^-1*b^-2*fresnels(a+b*x)+-1*pi^-1*b^-2*(b*c+-1*a*d)*sin((1/2)*pi*(a+b*x)^2)+-1//2*b^-2*d^-1*(b*c+-1*a*d)^2*fresnelc(a+b*x)+-1//2*d*pi^-1*b^-2*(a+b*x)*sin((1/2)*pi*(a+b*x)^2))
@test integrate(fresnelc(a+b*x), x) == :(b^-1*(a+b*x)*fresnelc(a+b*x)+-1*pi^-1*b^-1*sin((1/2)*pi*(a+b*x)^2))
@test integrate(x^3*fresnelc(a+b*x), x) == :(1//4*x^4*fresnelc(a+b*x)+-1//4*a^4*b^-4*fresnelc(a+b*x)+3//4*pi^-2*b^-4*fresnelc(a+b*x)+pi^-1*a^3*b^-4*sin((1/2)*pi*(a+b*x)^2)+2*a*pi^-2*b^-4*cos((1/2)*pi*(a+b*x)^2)+-1//4*pi^-1*b^-4*(a+b*x)^3*sin((1/2)*pi*(a+b*x)^2)+-1//4*pi^-2*b^-4*(3a+3*b*x)*cos((1/2)*pi*(a+b*x)^2)+3//2*pi^-1*a^2*b^-4*fresnels(a+b*x)+a*pi^-1*b^-4*(a+b*x)^2*sin((1/2)*pi*(a+b*x)^2)+-3//2*pi^-1*a^2*b^-4*(a+b*x)*sin((1/2)*pi*(a+b*x)^2))
@test integrate(x^2*fresnelc(a+b*x), x) == :(1//3*x^3*fresnelc(a+b*x)+-2//3*pi^-2*b^-3*cos((1/2)*pi*(a+b*x)^2)+1//3*a^3*b^-3*fresnelc(a+b*x)+-1*a*pi^-1*b^-3*fresnels(a+b*x)+-1*pi^-1*a^2*b^-3*sin((1/2)*pi*(a+b*x)^2)+-1//3*pi^-1*b^-3*(a+b*x)^2*sin((1/2)*pi*(a+b*x)^2)+a*pi^-1*b^-3*(a+b*x)*sin((1/2)*pi*(a+b*x)^2))
@test integrate(x*fresnelc(a+b*x), x) == :((1/2)*x^2*fresnelc(a+b*x)+(1/2)*pi^-1*b^-2*fresnels(a+b*x)+-1//2*a^2*b^-2*fresnelc(a+b*x)+a*pi^-1*b^-2*sin((1/2)*pi*(a+b*x)^2)+-1//2*pi^-1*b^-2*(a+b*x)*sin((1/2)*pi*(a+b*x)^2))
@test integrate(fresnelc(a+b*x), x) == :(b^-1*(a+b*x)*fresnelc(a+b*x)+-1*pi^-1*b^-1*sin((1/2)*pi*(a+b*x)^2))
@test integrate(x^7*fresnelc(b*x)^2, x) == :(1//8*x^8*fresnelc(b*x)^2+-10*pi^-5*b^-8*sin(pi*b^2*x^2)+-105//8*pi^-4*b^-8*fresnelc(b*x)^2+-105//16*pi^-4*b^-6*x^2+7//48*pi^-2*b^-2*x^6+-1//16*pi^-2*b^-2*x^6*cos(pi*b^2*x^2)+5//8*pi^-3*b^-4*x^4*sin(pi*b^2*x^2)+55//16*pi^-4*b^-6*x^2*cos(pi*b^2*x^2)+-7//4*pi^-2*b^-3*x^5*cos((1/2)*pi*b^2*x^2)*fresnelc(b*x)+-1//4*pi^-1*b^-1*x^7*fresnelc(b*x)*sin((1/2)*pi*b^2*x^2)+35//4*pi^-3*b^-5*x^3*fresnelc(b*x)*sin((1/2)*pi*b^2*x^2)+105//4*x*pi^-4*b^-7*cos((1/2)*pi*b^2*x^2)*fresnelc(b*x))
@test integrate(x^6*fresnelc(b*x)^2, x) == :(1//7*x^7*fresnelc(b*x)^2+-48//7*x*pi^-4*b^-6+6//35*pi^-2*b^-2*x^5+-531//112*2^(1/2)*pi^-4*b^-7*fresnelc(b*x*2^(1/2))+-1//14*pi^-2*b^-2*x^5*cos(pi*b^2*x^2)+17//28*pi^-3*b^-4*x^3*sin(pi*b^2*x^2)+21//8*x*pi^-4*b^-6*cos(pi*b^2*x^2)+96//7*pi^-4*b^-7*cos((1/2)*pi*b^2*x^2)*fresnelc(b*x)+-12//7*pi^-2*b^-3*x^4*cos((1/2)*pi*b^2*x^2)*fresnelc(b*x)+-2//7*pi^-1*b^-1*x^6*fresnelc(b*x)*sin((1/2)*pi*b^2*x^2)+48//7*pi^-3*b^-5*x^2*fresnelc(b*x)*sin((1/2)*pi*b^2*x^2))
@test integrate(x^4*fresnelc(b*x)^2, x) == :(1//5*x^5*fresnelc(b*x)^2+4//15*pi^-2*b^-2*x^3+-43//40*2^(1/2)*pi^-3*b^-5*fresnels(b*x*2^(1/2))+-1//10*pi^-2*b^-2*x^3*cos(pi*b^2*x^2)+11//20*x*pi^-3*b^-4*sin(pi*b^2*x^2)+16//5*pi^-3*b^-5*fresnelc(b*x)*sin((1/2)*pi*b^2*x^2)+-8//5*pi^-2*b^-3*x^2*cos((1/2)*pi*b^2*x^2)*fresnelc(b*x)+-2//5*pi^-1*b^-1*x^4*fresnelc(b*x)*sin((1/2)*pi*b^2*x^2))
@test integrate(x^3*fresnelc(b*x)^2, x) == :(1//4*x^4*fresnelc(b*x)^2+(1/2)*pi^-3*b^-4*sin(pi*b^2*x^2)+3//4*pi^-2*b^-4*fresnelc(b*x)^2+3//8*pi^-2*b^-2*x^2+-1//8*pi^-2*b^-2*x^2*cos(pi*b^2*x^2)+-3//2*x*pi^-2*b^-3*cos((1/2)*pi*b^2*x^2)*fresnelc(b*x)+-1//2*pi^-1*b^-1*x^3*fresnelc(b*x)*sin((1/2)*pi*b^2*x^2))
@test integrate(x^2*fresnelc(b*x)^2, x) == :(1//3*x^3*fresnelc(b*x)^2+2//3*x*pi^-2*b^-2+-4//3*pi^-2*b^-3*cos((1/2)*pi*b^2*x^2)*fresnelc(b*x)+-1//6*x*pi^-2*b^-2*cos(pi*b^2*x^2)+5//12*2^(1/2)*pi^-2*b^-3*fresnelc(b*x*2^(1/2))+-2//3*pi^-1*b^-1*x^2*fresnelc(b*x)*sin((1/2)*pi*b^2*x^2))
@test integrate(fresnelc(b*x)^2, x) == :(x*fresnelc(b*x)^2+(1/2)*2^(1/2)*pi^-1*b^-1*fresnels(b*x*2^(1/2))+-2*pi^-1*b^-1*fresnelc(b*x)*sin((1/2)*pi*b^2*x^2))
@test integrate(x^-5*fresnelc(b*x)^2, x) == :(-1//4*x^-4*fresnelc(b*x)^2+-1//24*b^2*x^-2+-1//12*pi*b^4*sinintegral(pi*b^2*x^2)+-1//12*pi^2*b^4*fresnelc(b*x)^2+-1//24*b^2*x^-2*cos(pi*b^2*x^2)+-1//6*b*x^-3*cos((1/2)*pi*b^2*x^2)*fresnelc(b*x)+1//6*pi*b^3*x^-1*fresnelc(b*x)*sin((1/2)*pi*b^2*x^2))
@test integrate(x^-9*fresnelc(b*x)^2, x) == :(-1//8*x^-8*fresnelc(b*x)^2+-1//336*b^2*x^-6+-1//336*b^2*x^-6*cos(pi*b^2*x^2)+1//280*pi^3*b^8*sinintegral(pi*b^2*x^2)+1//840*pi^4*b^8*fresnelc(b*x)^2+1//1680*pi^2*b^6*x^-2+-1//28*b*x^-7*cos((1/2)*pi*b^2*x^2)*fresnelc(b*x)+1//336*pi^2*b^6*x^-2*cos(pi*b^2*x^2)+1//420*pi*b^4*x^-4*sin(pi*b^2*x^2)+-1//420*pi^3*b^7*x^-1*fresnelc(b*x)*sin((1/2)*pi*b^2*x^2)+1//140*pi*b^3*x^-5*fresnelc(b*x)*sin((1/2)*pi*b^2*x^2)+1//420*pi^2*b^5*x^-3*cos((1/2)*pi*b^2*x^2)*fresnelc(b*x))
@test integrate(fresnelc(a+b*x)^2, x) == :(b^-1*fresnelc(a+b*x)^2*(a+b*x)+(1/2)*2^(1/2)*pi^-1*b^-1*fresnels(2^(1/2)*(a+b*x))+-2*pi^-1*b^-1*fresnelc(a+b*x)*sin((1/2)*pi*(a+b*x)^2))
@test integrate(x^2*fresnelc(d*(a+b*log(c*x^n))), x) == :(1//3*x^3*fresnelc(d*(a+b*log(c*x^n)))+x^3*(c*x^n)^(-3*n^-1)*(1//12+1//12*im)*erf(pi^-1//2*b^-1*d^-1*(1/2+(1/2)*im)*(3*n^-1+pi*im*a*b*d^2+pi*im*b^2*d^2*log(c*x^n)))*exp(-3*a*b^-1*n^-1+9//2*im*pi^-1*b^-2*d^-2*n^-2)+-1*x^3*(c*x^n)^(-3*n^-1)*(1//12+1//12*im)*erfi(pi^-1//2*b^-1*d^-1*(1/2+(1/2)*im)*(3*n^-1+-1*pi*im*a*b*d^2+-1*pi*im*b^2*d^2*log(c*x^n)))*exp(-3*a*b^-1*n^-1+-9//2*im*pi^-1*b^-2*d^-2*n^-2))
@test integrate(x*fresnelc(d*(a+b*log(c*x^n))), x) == :((1/2)*x^2*fresnelc(d*(a+b*log(c*x^n)))+x^2*(c*x^n)^(-2*n^-1)*(1//8+1//8*im)*erf(pi^-1//2*b^-1*d^-1*(1/2+(1/2)*im)*(2*n^-1+pi*im*a*b*d^2+pi*im*b^2*d^2*log(c*x^n)))*exp(pi^-1*b^-2*d^-2*n^-2*(2im+-2*pi*a*b*n*d^2))+-1*x^2*(c*x^n)^(-2*n^-1)*(1//8+1//8*im)*erfi(pi^-1//2*b^-1*d^-1*(1/2+(1/2)*im)*(2*n^-1+-1*pi*im*a*b*d^2+-1*pi*im*b^2*d^2*log(c*x^n)))*exp(-1*pi^-1*b^-2*d^-2*n^-2*(2im+2*pi*a*b*n*d^2)))
@test integrate(fresnelc(d*(a+b*log(c*x^n))), x) == :(x*fresnelc(d*(a+b*log(c*x^n)))+x*(c*x^n)^(-1*n^-1)*(1//4+1//4*im)*erf(pi^-1//2*b^-1*d^-1*(1/2+(1/2)*im)*(n^-1+pi*im*a*b*d^2+pi*im*b^2*d^2*log(c*x^n)))*exp(-1//2*b^-2*n^-2*(-1*im*pi^-1*d^-2+2*a*b*n))+-1*x*(c*x^n)^(-1*n^-1)*(1//4+1//4*im)*erfi(pi^-1//2*b^-1*d^-1*(1/2+(1/2)*im)*(n^-1+-1*pi*im*a*b*d^2+-1*pi*im*b^2*d^2*log(c*x^n)))*exp(-1//2*b^-2*n^-2*(im*pi^-1*d^-2+2*a*b*n)))
@test integrate(x^-1*fresnelc(d*(a+b*log(c*x^n))), x) == :(b^-1*n^-1*(a+b*log(c*x^n))*fresnelc(d*(a+b*log(c*x^n)))+-1*pi^-1*b^-1*d^-1*n^-1*sin((1/2)*pi*d^2*(a+b*log(c*x^n))^2))
@test integrate(x^-2*fresnelc(d*(a+b*log(c*x^n))), x) == :(-1*x^-1*fresnelc(d*(a+b*log(c*x^n)))+x^-1*(c*x^n)^(n^-1)*(1//4+1//4*im)*erf(pi^-1//2*b^-1*d^-1*(1/2+(1/2)*im)*(n^-1+-1*pi*im*a*b*d^2+-1*pi*im*b^2*d^2*log(c*x^n)))*exp((1/2)*b^-2*n^-2*(im*pi^-1*d^-2+2*a*b*n))+-1*x^-1*(c*x^n)^(n^-1)*(1//4+1//4*im)*erfi(pi^-1//2*b^-1*d^-1*(1/2+(1/2)*im)*(n^-1+pi*im*a*b*d^2+pi*im*b^2*d^2*log(c*x^n)))*exp((1/2)*b^-2*n^-2*(-1*im*pi^-1*d^-2+2*a*b*n)))
@test integrate(x^-3*fresnelc(d*(a+b*log(c*x^n))), x) == :(-1//2*x^-2*fresnelc(d*(a+b*log(c*x^n)))+x^-2*(c*x^n)^(2*n^-1)*(1//8+1//8*im)*erf(pi^-1//2*b^-1*d^-1*(1/2+(1/2)*im)*(2*n^-1+-1*pi*im*a*b*d^2+-1*pi*im*b^2*d^2*log(c*x^n)))*exp(pi^-1*b^-2*d^-2*n^-2*(2im+2*pi*a*b*n*d^2))+-1*x^-2*(c*x^n)^(2*n^-1)*(1//8+1//8*im)*erfi(pi^-1//2*b^-1*d^-1*(1/2+(1/2)*im)*(2*n^-1+pi*im*a*b*d^2+pi*im*b^2*d^2*log(c*x^n)))*exp(-1*pi^-1*b^-2*d^-2*n^-2*(2im+-2*pi*a*b*n*d^2)))
@test integrate((e*x)^m*fresnelc(d*(a+b*log(c*x^n))), x) == :(e^-1*(e*x)^(1+m)*(1+m)^-1*fresnelc(d*(a+b*log(c*x^n)))+x*(c*x^n)^(-1*n^-1*(1+m))*(e*x)^m*(1+m)^-1*(1//4+1//4*im)*erf(pi^-1//2*b^-1*d^-1*n^-1*(1/2+(1/2)*im)*(1+m+pi*im*a*b*n*d^2+pi*im*n*b^2*d^2*log(c*x^n)))*exp((1/2)*im*pi^-1*b^-2*d^-2*n^-2*(1+m)*(1+m+2*pi*im*a*b*n*d^2))+-1*x*(c*x^n)^(-1*n^-1*(1+m))*(e*x)^m*(1+m)^-1*(1//4+1//4*im)*erfi(pi^-1//2*b^-1*d^-1*n^-1*(1/2+(1/2)*im)*(1+m+-1*pi*im*a*b*n*d^2+-1*pi*im*n*b^2*d^2*log(c*x^n)))*exp(-1//2*im*pi^-1*b^-2*d^-2*n^-2*(1+m)*(1+m+-2*pi*im*a*b*n*d^2)))
@test integrate(fresnelc(b*x)^2*cos((1/2)*pi*b^2*x^2), x) == :(1//3*b^-1*fresnelc(b*x)^3)
@test integrate(cos((1/2)*pi*b^2*x^2)*fresnelc(b*x), x) == :((1/2)*b^-1*fresnelc(b*x)^2)
@test integrate(fresnelc(b*x)^-1*cos((1/2)*pi*b^2*x^2), x) == :(b^-1*log(fresnelc(b*x)))
@test integrate(fresnelc(b*x)^-2*cos((1/2)*pi*b^2*x^2), x) == :(-1*b^-1*fresnelc(b*x)^-1)
@test integrate(fresnelc(b*x)^-3*cos((1/2)*pi*b^2*x^2), x) == :(-1//2*b^-1*fresnelc(b*x)^-2)
@test integrate(fresnelc(b*x)^n*cos((1/2)*pi*b^2*x^2), x) == :(b^-1*(1+n)^-1*fresnelc(b*x)^(1+n))
@test integrate(x^8*cos((1/2)*pi*b^2*x^2)*fresnelc(b*x), x) == :(40*pi^-5*b^-9*sin(pi*b^2*x^2)+-7//12*pi^-2*b^-3*x^6+105//2*pi^-4*b^-9*fresnelc(b*x)^2+105//4*pi^-4*b^-7*x^2+-55//4*pi^-4*b^-7*x^2*cos(pi*b^2*x^2)+-5//2*pi^-3*b^-5*x^4*sin(pi*b^2*x^2)+1//4*pi^-2*b^-3*x^6*cos(pi*b^2*x^2)+pi^-1*b^-2*x^7*fresnelc(b*x)*sin((1/2)*pi*b^2*x^2)+-105*x*pi^-4*b^-8*cos((1/2)*pi*b^2*x^2)*fresnelc(b*x)+-35*pi^-3*b^-6*x^3*fresnelc(b*x)*sin((1/2)*pi*b^2*x^2)+7*pi^-2*b^-4*x^5*cos((1/2)*pi*b^2*x^2)*fresnelc(b*x))
@test integrate(x^7*cos((1/2)*pi*b^2*x^2)*fresnelc(b*x), x) == :(24*x*pi^-4*b^-7+-3//5*pi^-2*b^-3*x^5+-48*pi^-4*b^-8*cos((1/2)*pi*b^2*x^2)*fresnelc(b*x)+-147//16*x*pi^-4*b^-7*cos(pi*b^2*x^2)+-17//8*pi^-3*b^-5*x^3*sin(pi*b^2*x^2)+1//4*pi^-2*b^-3*x^5*cos(pi*b^2*x^2)+531//32*2^(1/2)*pi^-4*b^-8*fresnelc(b*x*2^(1/2))+pi^-1*b^-2*x^6*fresnelc(b*x)*sin((1/2)*pi*b^2*x^2)+-24*pi^-3*b^-6*x^2*fresnelc(b*x)*sin((1/2)*pi*b^2*x^2)+6*pi^-2*b^-4*x^4*cos((1/2)*pi*b^2*x^2)*fresnelc(b*x))
@test integrate(x^5*cos((1/2)*pi*b^2*x^2)*fresnelc(b*x), x) == :(-2//3*pi^-2*b^-3*x^3+-8*pi^-3*b^-6*fresnelc(b*x)*sin((1/2)*pi*b^2*x^2)+-11//8*x*pi^-3*b^-5*sin(pi*b^2*x^2)+1//4*pi^-2*b^-3*x^3*cos(pi*b^2*x^2)+43//16*2^(1/2)*pi^-3*b^-6*fresnels(b*x*2^(1/2))+pi^-1*b^-2*x^4*fresnelc(b*x)*sin((1/2)*pi*b^2*x^2)+4*pi^-2*b^-4*x^2*cos((1/2)*pi*b^2*x^2)*fresnelc(b*x))
@test integrate(x^4*cos((1/2)*pi*b^2*x^2)*fresnelc(b*x), x) == :(-1*pi^-3*b^-5*sin(pi*b^2*x^2)+-3//2*pi^-2*b^-5*fresnelc(b*x)^2+-3//4*pi^-2*b^-3*x^2+1//4*pi^-2*b^-3*x^2*cos(pi*b^2*x^2)+pi^-1*b^-2*x^3*fresnelc(b*x)*sin((1/2)*pi*b^2*x^2)+3*x*pi^-2*b^-4*cos((1/2)*pi*b^2*x^2)*fresnelc(b*x))
@test integrate(x^3*cos((1/2)*pi*b^2*x^2)*fresnelc(b*x), x) == :(-1*x*pi^-2*b^-3+2*pi^-2*b^-4*cos((1/2)*pi*b^2*x^2)*fresnelc(b*x)+-5//8*2^(1/2)*pi^-2*b^-4*fresnelc(b*x*2^(1/2))+1//4*x*pi^-2*b^-3*cos(pi*b^2*x^2)+pi^-1*b^-2*x^2*fresnelc(b*x)*sin((1/2)*pi*b^2*x^2))
@test integrate(x*cos((1/2)*pi*b^2*x^2)*fresnelc(b*x), x) == :(pi^-1*b^-2*fresnelc(b*x)*sin((1/2)*pi*b^2*x^2)+-1//4*2^(1/2)*pi^-1*b^-2*fresnels(b*x*2^(1/2)))
@test integrate(cos((1/2)*pi*b^2*x^2)*fresnelc(b*x), x) == :((1/2)*b^-1*fresnelc(b*x)^2)
@test integrate(x^-4*cos((1/2)*pi*b^2*x^2)*fresnelc(b*x), x) == :(-1//12*b*x^-2+-1//3*x^-3*cos((1/2)*pi*b^2*x^2)*fresnelc(b*x)+-1//6*pi*b^3*sinintegral(pi*b^2*x^2)+-1//6*pi^2*b^3*fresnelc(b*x)^2+-1//12*b*x^-2*cos(pi*b^2*x^2)+1//3*pi*b^2*x^-1*fresnelc(b*x)*sin((1/2)*pi*b^2*x^2))
@test integrate(x^-8*cos((1/2)*pi*b^2*x^2)*fresnelc(b*x), x) == :(-1//84*b*x^-6+-1//7*x^-7*cos((1/2)*pi*b^2*x^2)*fresnelc(b*x)+-1//84*b*x^-6*cos(pi*b^2*x^2)+1//70*pi^3*b^7*sinintegral(pi*b^2*x^2)+1//210*pi^4*b^7*fresnelc(b*x)^2+1//420*pi^2*b^5*x^-2+1//84*pi^2*b^5*x^-2*cos(pi*b^2*x^2)+1//105*pi*b^3*x^-4*sin(pi*b^2*x^2)+-1//105*pi^3*b^6*x^-1*fresnelc(b*x)*sin((1/2)*pi*b^2*x^2)+1//35*pi*b^2*x^-5*fresnelc(b*x)*sin((1/2)*pi*b^2*x^2)+1//105*pi^2*b^4*x^-3*cos((1/2)*pi*b^2*x^2)*fresnelc(b*x))
@test integrate(x^7*fresnelc(b*x)*sin((1/2)*pi*b^2*x^2), x) == :(-4*pi^-3*b^-5*x^3+1//14*pi^-1*b^-1*x^7+-48*pi^-4*b^-8*fresnelc(b*x)*sin((1/2)*pi*b^2*x^2)+-147//16*x*pi^-4*b^-7*sin(pi*b^2*x^2)+1//4*pi^-2*b^-3*x^5*sin(pi*b^2*x^2)+17//8*pi^-3*b^-5*x^3*cos(pi*b^2*x^2)+531//32*2^(1/2)*pi^-4*b^-8*fresnels(b*x*2^(1/2))+-1*pi^-1*b^-2*x^6*cos((1/2)*pi*b^2*x^2)*fresnelc(b*x)+6*pi^-2*b^-4*x^4*fresnelc(b*x)*sin((1/2)*pi*b^2*x^2)+24*pi^-3*b^-6*x^2*cos((1/2)*pi*b^2*x^2)*fresnelc(b*x))
@test integrate(x^6*fresnelc(b*x)*sin((1/2)*pi*b^2*x^2), x) == :(-15//2*pi^-3*b^-7*fresnelc(b*x)^2+-15//4*pi^-3*b^-5*x^2+-11//2*pi^-4*b^-7*sin(pi*b^2*x^2)+1//12*pi^-1*b^-1*x^6+1//4*pi^-2*b^-3*x^4*sin(pi*b^2*x^2)+7//4*pi^-3*b^-5*x^2*cos(pi*b^2*x^2)+-1*pi^-1*b^-2*x^5*cos((1/2)*pi*b^2*x^2)*fresnelc(b*x)+5*pi^-2*b^-4*x^3*fresnelc(b*x)*sin((1/2)*pi*b^2*x^2)+15*x*pi^-3*b^-6*cos((1/2)*pi*b^2*x^2)*fresnelc(b*x))
@test integrate(x^5*fresnelc(b*x)*sin((1/2)*pi*b^2*x^2), x) == :(-4*x*pi^-3*b^-5+1//10*pi^-1*b^-1*x^5+8*pi^-3*b^-6*cos((1/2)*pi*b^2*x^2)*fresnelc(b*x)+-43//16*2^(1/2)*pi^-3*b^-6*fresnelc(b*x*2^(1/2))+1//4*pi^-2*b^-3*x^3*sin(pi*b^2*x^2)+11//8*x*pi^-3*b^-5*cos(pi*b^2*x^2)+-1*pi^-1*b^-2*x^4*cos((1/2)*pi*b^2*x^2)*fresnelc(b*x)+4*pi^-2*b^-4*x^2*fresnelc(b*x)*sin((1/2)*pi*b^2*x^2))
@test integrate(x^3*fresnelc(b*x)*sin((1/2)*pi*b^2*x^2), x) == :(1//6*pi^-1*b^-1*x^3+2*pi^-2*b^-4*fresnelc(b*x)*sin((1/2)*pi*b^2*x^2)+-5//8*2^(1/2)*pi^-2*b^-4*fresnels(b*x*2^(1/2))+1//4*x*pi^-2*b^-3*sin(pi*b^2*x^2)+-1*pi^-1*b^-2*x^2*cos((1/2)*pi*b^2*x^2)*fresnelc(b*x))
@test integrate(x^2*fresnelc(b*x)*sin((1/2)*pi*b^2*x^2), x) == :((1/2)*pi^-1*b^-3*fresnelc(b*x)^2+1//4*pi^-1*b^-1*x^2+1//4*pi^-2*b^-3*sin(pi*b^2*x^2)+-1*x*pi^-1*b^-2*cos((1/2)*pi*b^2*x^2)*fresnelc(b*x))
@test integrate(x*fresnelc(b*x)*sin((1/2)*pi*b^2*x^2), x) == :((1/2)*x*pi^-1*b^-1+-1*pi^-1*b^-2*cos((1/2)*pi*b^2*x^2)*fresnelc(b*x)+1//4*2^(1/2)*pi^-1*b^-2*fresnelc(b*x*2^(1/2)))
@test integrate(x^-2*fresnelc(b*x)*sin((1/2)*pi*b^2*x^2), x) == :(1//4*b*sinintegral(pi*b^2*x^2)+(1/2)*pi*b*fresnelc(b*x)^2+-1*x^-1*fresnelc(b*x)*sin((1/2)*pi*b^2*x^2))
@test integrate(x^-6*fresnelc(b*x)*sin((1/2)*pi*b^2*x^2), x) == :(-7//120*pi^2*b^5*sinintegral(pi*b^2*x^2)+-1//5*x^-5*fresnelc(b*x)*sin((1/2)*pi*b^2*x^2)+-1//30*pi^3*b^5*fresnelc(b*x)^2+-1//40*b*x^-4*sin(pi*b^2*x^2)+-1//60*pi*b^3*x^-2+-1//24*pi*b^3*x^-2*cos(pi*b^2*x^2)+-1//15*pi*b^2*x^-3*cos((1/2)*pi*b^2*x^2)*fresnelc(b*x)+1//15*pi^2*b^4*x^-1*fresnelc(b*x)*sin((1/2)*pi*b^2*x^2))
@test integrate(x^-10*fresnelc(b*x)*sin((1/2)*pi*b^2*x^2), x) == :(-1//9*x^-9*fresnelc(b*x)*sin((1/2)*pi*b^2*x^2)+-1//144*b*x^-8*sin(pi*b^2*x^2)+-1//756*pi*b^3*x^-6+1//1890*pi^5*b^9*fresnelc(b*x)^2+1//3780*pi^3*b^7*x^-2+83//30240*pi^4*b^9*sinintegral(pi*b^2*x^2)+-11//3024*pi*b^3*x^-6*cos(pi*b^2*x^2)+5//2016*pi^3*b^7*x^-2*cos(pi*b^2*x^2)+67//30240*pi^2*b^5*x^-4*sin(pi*b^2*x^2)+-1//63*pi*b^2*x^-7*cos((1/2)*pi*b^2*x^2)*fresnelc(b*x)+-1//945*pi^4*b^8*x^-1*fresnelc(b*x)*sin((1/2)*pi*b^2*x^2)+1//315*pi^2*b^4*x^-5*fresnelc(b*x)*sin((1/2)*pi*b^2*x^2)+1//945*pi^3*b^6*x^-3*cos((1/2)*pi*b^2*x^2)*fresnelc(b*x))

max_stars_repo_path: test/8 Special functions/8.2 Fresnel integral functions.jl


text:
mutable struct EVinstance
    arrivalTime::Float64
    departureTime::Float64
    reportedDepartureTime::Float64   ##para el caso en que hay incertidumbre
    requestedEnergy::Float64
    chargingPower::Float64
    currentWorkload::Float64
    currentDeadline::Float64
    currentReportedDeadline::Float64
    currentPower::Float64
    departureWorkload::Float64
    completionTime::Float64

    #inicializo la instancia solo con los 4 importantes y completo los otros al comienzo.
    #En particular reportedDepartureTime lo pongo igual a Departure time por defecto.
    #Departure workload y CompletionTime queda en NaN hasta que se calculen mas adelante.
    EVinstance( arrivalTime::Float64,
                departureTime::Float64,
                requestedEnergy::Float64,
                chargingPower::Float64) = new(  arrivalTime,
                                                departureTime,
                                                departureTime,
                                                requestedEnergy,
                                                chargingPower,
                                                requestedEnergy,
                                                departureTime-arrivalTime,
                                                departureTime-arrivalTime,
                                                0.0,NaN,NaN)

    #este segundo constructor me deja llenar el reportedDepartureTime
    #en particular el valor de currentReportedDeadline se calcula con lo el reportado
    EVinstance( arrivalTime::Float64,
                departureTime::Float64,
                reportedDepartureTime::Float64,
                requestedEnergy::Float64,
                chargingPower::Float64) = new(  arrivalTime,
                                                departureTime,
                                                reportedDepartureTime,
                                                requestedEnergy,
                                                chargingPower,
                                                requestedEnergy,
                                                departureTime-arrivalTime,
                                                reportedDepartureTime-arrivalTime,
                                                0.0,NaN,NaN)
end

mutable struct Snapshot
    t::Float64                          #snapshot time
    charging::Array{EVinstance}         #vehicles in the system in charge
    alreadyCharged::Array{EVinstance}   #already charged vehicles still present
end

mutable struct TimeTrace
    T::Vector{Float64}          #event times
    X::Vector{UInt16}           #charging vehicles
    Y::Vector{UInt16}           #already charged
    P::Vector{Float64}          #used power
end

mutable struct SimStatistics
    rangeX::Vector{Integer}
    pX::Vector{Float64} #steady state X
    rangeY::Vector{Integer}
    pY::Vector{Float64} #steady state Y
    avgX::Float64       #average X
    avgY::Float64       #average Y
    pD::Float64         #probability of expired deadline
    avgW::Float64       #average unfinished workload (taking finished into account)
end

#defino la estructura resultados de simulacion
mutable struct EVSim
    parameters::Dict
    timetrace::TimeTrace
    EVs::Vector{EVinstance}
    snapshots::Vector{Snapshot}
    stats::SimStatistics
end


text:
function solve()
    x, y, a, b = [parse(BigInt, x) for x in split(readline())]
    ret = 0
    while a*x ≤ x+b && a*x<y
        x *= a
        ret += 1
    end
    ret + (y-1-x)÷b
end

println(solve())

max_stars_repo_name: aishikawa/atcoder-julia

max_issues_repo_name: aishikawa/atcoder-julia

max_forks_repo_name: aishikawa/atcoder-julia


text:
szeros(::Type{T}, N) where T = @SVector zeros(T, N)
szeros(N)= @SVector zeros(N)
szeros(::Type{T}, N1, N2) where T = @SMatrix zeros(T, N1, N2)
szeros(N1, N2)= @SMatrix zeros(N1, N2)

sones(::Type{T}, N) where T = @SVector ones(T, N)
sones(N)= @SVector ones(N)
sones(::Type{T}, N1, N2) where T = @SMatrix ones(T, N1, N2)
sones(N1, N2)= @SMatrix ones(N1, N2)

srand(::Type{T}, N) where T = @SVector rand(T, N)
srand(N)= @SVector rand(N)
srand(::Type{T}, N1, N2) where T = @SMatrix rand(T, N1, N2)
srand(N1, N2)= @SMatrix rand(N1, N2)

# To fix StaticArray bug
Base.:*(::SMatrix{N,0,T,0},::SVector{0}) where {T,N} = szeros(T,N)
Base.:*(::SMatrix{N,0,T,0},::SMatrix{0,M,T,0}) where {T,N,M} = szeros(T,N,M)
Base.inv(::SMatrix{0,0,T,0}) where T = szeros(T,0,0)


text:
printstyled("  Utils\n", color=:light_green)
printstyled("    info_dump\n", color=:light_green)

redirect_stdout(out) do
  info_dump(U)
end


text:
@testset "RichardsEquation" begin

    @testset "HF135" begin

    ϵ = 1e-7

    # Wetting of an HF135 membrane, Van Genuchten model
    # Data from Buser (PhD thesis, 2016)
    # http://hdl.handle.net/1773/38064
    θr = 0.0473
    θs = 0.945
    k = 5.50e-13  # m**2
    α = 0.2555  # 1/m
    n = 2.3521
    θi = 0.102755  # Computed from P0

    θb = θs - ϵ

    model = VanGenuchten(n=n, α=α, k=k, θr=θr, θs=θs)

    dprob = DirichletProblem(DiffusionEquation(model), i=θi, b=θb)

    θ = solve(dprob, itol=1e-7)

    hb = 0
    hi = -30.591486389337

    prob = DirichletProblem(RichardsEquation(model), i=hi, b=hb)

    h = solve(prob)

    r = range(0, 0.05, length=500)
    t = [10 2.73 314]

    @test all(@. isapprox(θ(r,t), θh(model, h(r,t)), atol=1e-1))
    @test all(@. isapprox(flux(θ,r[2:end],t), flux(h,r[2:end],t), atol=1e-4))
    end

end


text:
# Autogenerated wrapper script for LLVM_assert_jll for aarch64-linux-musl-cxx11
export liblto, llc, lld, llvm_mca, opt

using Zlib_jll
using libLLVM_assert_jll
JLLWrappers.@generate_wrapper_header("LLVM_assert")
JLLWrappers.@declare_library_product(liblto, "libLTO.so.13jl")
JLLWrappers.@declare_executable_product(llc)
JLLWrappers.@declare_executable_product(lld)
JLLWrappers.@declare_executable_product(llvm_mca)
JLLWrappers.@declare_executable_product(opt)
function __init__()
    JLLWrappers.@generate_init_header(Zlib_jll, libLLVM_assert_jll)
    JLLWrappers.@init_library_product(
        liblto,
        "lib/libLTO.so",
        nothing,
    )

    JLLWrappers.@init_executable_product(
        llc,
        "tools/llc",
    )

    JLLWrappers.@init_executable_product(
        lld,
        "bin/lld",
    )

    JLLWrappers.@init_executable_product(
        llvm_mca,
        "tools/llvm-mca",
    )

    JLLWrappers.@init_executable_product(
        opt,
        "tools/opt",
    )

    JLLWrappers.@generate_init_footer()
end  # __init__()


text:
using HypothesisTests, Test
using HypothesisTests: default_tail

@testset "F-tests" begin
    @testset "Basic variance F-test" begin
        Random.seed!(12)
        y1_h0 = 4 .+ randn(500)
        y2_h0 = 4 .+ randn(400)

        t = VarianceFTest(y1_h0, y2_h0)

        @test t.n_x == 500
        @test t.n_y == 400
        @test t.df_x == 499
        @test t.df_y == 399
        @test t.F ≈ 0.974693 rtol=1e-5
        @test pvalue(t) ≈ 0.784563 rtol=1e-5
        @test pvalue(t, tail=:left) ≈ 0.392281 rtol=1e-5
        @test pvalue(t, tail=:right) ≈ 0.607718 rtol=1e-5
        @test default_tail(t) == :both

        t = VarianceFTest(y2_h0, y1_h0)

        @test t.n_x == 400
        @test t.n_y == 500
        @test t.df_x == 399
        @test t.df_y == 499
        @test t.F ≈ 1.025963 rtol=1e-5
        @test pvalue(t) ≈ 0.784563 rtol=1e-5
        @test pvalue(t, tail=:right) ≈ 0.392281 rtol=1e-5
        @test pvalue(t, tail=:left) ≈ 0.607718 rtol=1e-5
        @test default_tail(t) == :both

        y1_h1 = 0.7*randn(200)
        y2_h1 = 1.3*randn(120)

        t = VarianceFTest(y1_h1, y2_h1)

        @test t.n_x == 200
        @test t.n_y == 120
        @test t.df_x == 199
        @test t.df_y == 119
        @test t.F ≈ 0.367547 rtol=1e-5
        @test pvalue(t) < 1e-8
        @test default_tail(t) == :both
        @test pvalue(t, tail=:left) < 1e-8
        @test pvalue(t, tail=:right) > 1.0 - 1e-8
    end
end


text:
using LinearAlgebra
using SparseArrays

"""
Construct a Heisenberg-type Hamiltonian as sparse matrix using site-dependent interaction strength, with nearest-neighbor interactions according to adjacency.
"""
function construct_hamiltonian(J::Array{<:Real,3})
    L = size(J, 3)
    @assert(size(J) == (3, L, L))

    # spin operators (Pauli matrices divided by 2)
    sigma12 = (0.5*sparse([0.  1.; 1.  0.]),
               0.5*sparse([0. -im; im  0.]),
               0.5*sparse([1.  0.; 0. -1.]))

    H = spzeros(Float64, 2^L, 2^L)

    # interaction terms
    for i in 1:L
        for j in i+1:L
            for k in 1:3
                # considering only entries in J for i < j
                if J[k, i, j] ≠ 0
                    H -= J[k, i, j] * real(kron(sparse_identity(2^(L-j)), sigma12[k], sparse_identity(2^(j-i-1)), sigma12[k], sparse_identity(2^(i-1))))
                end
            end
        end
    end

    return H
end

max_stars_repo_name: keefehuang/vqt

max_issues_repo_name: keefehuang/vqt

max_forks_repo_name: keefehuang/vqt


text:
# addprocs, same fun everywhere!
using Distributed
while nprocs()<4
    addprocs(1)
end
@everywhere using Pkg
@everywhere Pkg.activate(".")
using DEMC

@everywhere using ParallelDataTransfer
@everywhere using DEMC
#@everywhere include("./src/DEMC.jl")
@everywhere using Random
@everywhere Random.seed!(31953150 + myid())
# set up target distribution: Multivariate Normal
ndim = 30 # Number of dimensions
μ = rand(ndim) # mean of each dimension
passobj(myid(), workers(), :μ)
# log objective function
@everywhere log_obj(mean) = -sum((mean.-μ).^2)

# set up of DEMCz chain
Npar = length(μ)
# set up of DEMCz chain
opts = DEMC.demcopt(Npar)
opts.blockindex = [1:Npar] # parameter blocks: here choose all parameters to be updated simultaenously
opts.Nblocks = length(opts.blockindex)
opts.eps_scale = 1e-5*ones(Npar) # scale of random error around DE update
opts.γ = 2.38 # scale of DE update, 2.38 is the "optimal" number for a normal distribution
opts.N = nworkers() # number of chains
opts.K = 10 # every K steps add current N draws to Z
# Number of iterations in Chain
opts.Ngeneration = 20_000
opts.TN = 0.0
opts.T0 = 5.0

Z = randn((10*Npar, Npar))*10 # initial distribution (completely off to make a difficult test case)

# Number of iterations in Chain
mc, Z = DEMC.demcz_anneal_par(log_obj, Z, opts; sync_every=500)

# check best element
bestval, bestel = findmax(mc.log_obj)
bestpar = mc.chain[bestel[1], :, bestel[2]]

@test abs(bestval)>-1e-1

#bestvals = [maximum(mc.log_obj[:, ig]) for ig = 1:Ngen]
#DEMC.plot(1:Ngen, bestvals)


text:
# BlockOperator supports adding rows and columns.

immutable BlockOperator{O,T} <: Operator{T}
    mat11::Matrix{T}
    mat12::Matrix{T}
    mat21::Matrix{T}
    op::O

    function BlockOperator(mat11::Matrix{T},mat12::Matrix{T},mat21::Matrix{T},op::O)
        @assert size(mat11,1)==size(mat12,1)
        @assert size(mat11,2)==size(mat21,2)
        new(mat11,mat12,mat21,op)
    end
end

BlockOperator(mat11::Matrix,mat12::Matrix,mat21::Matrix,B::Operator)=BlockOperator{typeof(B),
                                                               promote_type(eltype(mat11),eltype(mat12),eltype(mat21),
                                                                            eltype(B))}(mat11,mat12,mat21,B)

BlockOperator(mat11::Matrix,B::Operator)=BlockOperator(mat11,Array(eltype(mat11),size(mat11,1),0),
                                                            Array(eltype(mat11),0,size(mat11,2)),B)

function Base.hcat{S<:Number}(cols::Matrix{S},B::Operator)
    T = promote_type(S,eltype(B))
    BlockOperator{typeof(B),T}(Array(T,0,size(cols,2)),Array(T,0,0),cols,B)
end

Base.hcat{T<:Number}(cols::Vector{T},B::Operator)=hcat(reshape(cols,length(cols),1),B)

function BlockOperator{BO<:Operator}(A::Matrix{BO})
    @assert size(A,1)==1
    M=vec(A[1,1:end-1])
    B=A[1,end]

    rs=rangespace(B)

    if isa(rs,UnsetSpace)
        rs=UnsetSpace()  # for the case of all constants

        for k=1:size(A,2)-1
            if isa(A[1,k],Multiplication)
                # cols are going to be constantspace, so rangespace is just space of f
                rs=union(rs,space(A[1,k].f))
            end
        end

        if !isambiguous(rs)
            B=promotedomainspace(B,choosedomainspace(B,rs))
            rs=rangespace(B)
        else
            rs=UnsetSpace()
        end
    end

    T=mapreduce(eltype,promote_type,M)
    colsv=Array(Vector{T},length(M))
    for k=1:length(M)
        if isa(M[k],Multiplication)
            ds=domainspace(M[k])
            @assert isa(ds,UnsetSpace) || isa(ds,ConstantSpace)
            @assert !isambiguous(rs)
            colsv[k]=coefficients(M[k].f,rs)
        elseif isconstop(M[k])
            # if rs is UnsetSpace this just generates constantspace
            colsv[k]=Number(M[k])*ones(rs).coefficients
        else
            error("Not implemented")
        end
    end

    m=mapreduce(length,max,colsv)
    cols=zeros(T,m,length(M))
    for k=1:length(M)
        cols[1:length(colsv[k]),k]=colsv[k]
    end

    hcat(cols,B)
end

Base.convert{T}(::Type{Operator{T}},A::BlockOperator)=BlockOperator(convert(Matrix{T},A.mat11),
                        convert(Matrix{T},A.mat12),convert(Matrix{T},A.mat21),convert(Operator{T},A.op))

function rangespace(B::BlockOperator)
    rs=rangespace(B.op)
    if isa(rs,UnsetSpace) || size(B.mat11,1)==0
        rs # avoids TupleSpace⊕UnsetSpace
    else
        TupleSpace(fill(ConstantSpace(),size(B.mat11,1))...,rs)
    end
end

function domainspace(B::BlockOperator)
    ds=domainspace(B.op)
    if isa(ds,UnsetSpace) || size(B.mat11,2)==0
        ds # avoids TupleSpace⊕UnsetSpace
    else
        TupleSpace(fill(ConstantSpace(),size(B.mat11,2))...,ds)
    end
end

bandinds(B::BlockOperator)=min(1-size(B.mat21,1)-size(B.mat11,1),
                               bandinds(B.op,1)+size(B.mat11,2)-size(B.mat11,1)),
                           max(size(B.mat11,2)+size(B.mat12,2)-1,
                               bandinds(B.op,2)+size(B.mat11,2)-size(B.mat11,1))

function getindex(B::BlockOperator,k::Integer,j::Integer)
    n,m=size(B.mat11)
    if k ≤ n && j ≤ m
        B.mat11[k,j]
    elseif k ≤ n && j ≤ m + size(B.mat12,2)
        B.mat12[k,j-m]
    elseif k ≤ n + size(B.mat21,1) && j ≤ m
        B.mat21[k-n,j]
    elseif k > n && j > m
        B.op[k-n,j-m]
    else
        zero(eltype(B))
    end
end

choosedomainspace(B::BlockOperator,f::Space) =
    TupleSpace(fill(ConstantSpace(),size(B.mat11,2))...,choosedomainspace(B.op,f))

function promotedomainspace(P::BlockOperator,S::TupleSpace)
    m=size(P.mat11,2)
    @assert length(S.spaces)==m+1  #TODO: extra tuple?
    @assert all(sp->isa(sp,ConstantSpace),S.spaces[1:m])
    sp=S.spaces[end]

    op=promotedomainspace(P.op,sp)
    if size(P.mat11,1)==0 && size(P.mat21,1)==1 && isa(rangespace(P),UnsetSpace)
        # this is to allow unset space,
        # so we pass to the standard constructor
        BlockOperator([P.mat21 op])
    else
        # we don't know how to change the rangespace
        # TODO: convert coefficients from old rangespace to
        # new rangespace
        @assert rangespace(op)==rangespace(P.op)
        BlockOperator(P.mat11,P.mat12,P.mat21,op)
    end
end

# function promoterangespace(P::BlockOperator,S::TupleSpace)
#     m=size(P.mat11,1)
#     @assert length(S.spaces)==m+1  #TODO: extra tuple?
#     @assert all(sp->isa(sp,ConstantSpace),S.spaces[1:m])
#     sp=S.spaces[end]
#
#     op=promotedomainspace(P.op,sp)
#     if size(P.mat11,1)==0 && size(P.mat21,1)==1 && isa(rangespace(P),UnsetSpace)
#         # this is to allow unset space,
#         # so we pass to the standard constructor
#         BlockOperator([P.mat21 op])
#     else
#         # we don't know how to change the rangespace
#         # TODO: convert coefficients from old rangespace to
#         # new rangespace
#         @assert rangespace(op)==rangespace(P)
#         BlockOperator(P.mat11,P.mat12,P.mat21,op)
#     end
# end

## BlockFunctional

immutable BlockFunctional{T<:Number,B<:Operator} <: Operator{T}
    cols::Vector{T}
    op::B
end

@functional BlockFunctional

BlockFunctional{T<:Number}(cols::Vector{T},op::Operator) = BlockFunctional{promote_type(T,eltype(op)),typeof(op)}(promote_type(T,eltype(op))[cols],op)
BlockFunctional{T<:Number}(col::T,op::Operator) = BlockFunctional{promote_type(T,eltype(op)),typeof(op)}(promote_type(T,eltype(op))[col],op)

function domainspace(B::BlockFunctional)
    ds=domainspace(B.op)
    if isa(ds,UnsetSpace) || length(B.cols)==0
        ds # avoids TupleSpace⊕UnsetSpace
    else
        TupleSpace(fill(ConstantSpace(),length(B.cols))...,ds)
    end
end

function promotedomainspace(P::BlockFunctional,S::TupleSpace)
    @assert isa(S.spaces[1],ConstantSpace)
    sp=length(S.spaces)==2?S.spaces[2]:TupleSpace(S.spaces[2:end])

    op=promotedomainspace(P.op,sp)

    BlockFunctional(P.cols,op)
end

function Base.getindex{T<:Number}(P::BlockFunctional{T},kr::Range)
    lcols = length(P.cols)
    if kr[end]≤lcols
        P.cols[kr]
    else
        opr = intersect(kr,length(P.cols)+1:kr[end])
        [P.cols[intersect(kr,1:length(P.cols))];
         P.op[opr[1]-lcols:opr[end]-lcols]]
    end
end

function Base.getindex{T<:Number}(P::BlockFunctional{T},k::Integer)
    lcols = length(P.cols)
    if k≤lcols
        P.cols[k]
    else
        P.op[k-lcols]
    end
end

@eval Base.convert{T}(::Type{Operator{T}},P::BlockFunctional) =
    BlockFunctional(convert(Vector{T},P.cols),convert(Operator{T},P.op))


text:
print("Hello World! ☺ ♫ ☕ \n\n")

using Plots
gr()

plot([1,2],[2,3])

module rrt
    export Node, Edge, Obstacle, Room, Point

    struct Point
        x::Int
        y::Int
    end

    struct Node
        id::Int
        iPrev::Int #parent
        state::Point
    end

    struct Edge
        startnode::Int
        endnode::Int
    end

    struct Obstacle
        SW::Point
        NE::Point
    end

    ### Some types?
    struct Room #corners of the room
        x1::Int
        y1::Int #bottom left corner?
        x2::Int
        y2::Int
        obstacleList::Array{Obstacle}
    end

    struct robot
    end

end

function isCollidingNode(pt,obs)
    # this does not work.
    px,py = pt.x, pt.y
    x1,y1 = obs.SW.x, obs.SW.y
    x2,y2 = obs.NE.x, obs.NE.y

    if (px > x1 && px < x2 && py > y1 && py < y2)
        #print("Node in obstacle, discarded.")
        return true
    else
        return false
    end
end

function isCollidingEdge(r, nn, obs)
    # to detect collision, let's just check whether any of the four
    # lines of the rectangular obstacle intersect with our edge
    # ignore collinearity for now

    # Tofix: To make it look prettier in the graph, I am just going to add a 1 grid pt
    # spacer for now, until I figure out how to check for coincidence of point
    # in line. Or really, I should use GeometryShapes library

    x1,y1 = obs.SW.x, obs.SW.y
    x2,y2 = obs.NE.x, obs.NE.y

    x1 -= 1
    y1 -= 1
    x2 += 1
    y2 += 1

    pt1 = rrt.Point(x1, y1)
    pt2 = rrt.Point(x1, y2)
    pt3 = rrt.Point(x2, y2)
    pt4 = rrt.Point(x2, y1)

    # let A, B be r, nn
    coll1 = intersectLineSeg(r, nn, pt1, pt2)
    coll2 = intersectLineSeg(r, nn, pt2, pt3)
    coll3 = intersectLineSeg(r, nn, pt3, pt4)
    coll4 = intersectLineSeg(r, nn, pt4, pt1)

    if (!coll1 && !coll2 && !coll3 && !coll4)
        return false
    else
        #print("you're colliding!\n")
        #@show nn
        return true
    end
end

function inGoalRegion(node, goal)
    goal = rrt.Point(18,18)
    n = node.state
    if node == goal
        return true
    end
    distx = abs(n.x - goal.x)
    disty = abs(n.y - goal.y)
    if ( distx < 2 && disty < 2) #radius of goal
        return true
    end
    return false

end

function ccw(A,B,C)
    # determines direction of lines formed by three points
	return (C.y-A.y) * (B.x-A.x) > (B.y-A.y) * (C.x-A.x)
end

function intersectLineSeg(A,B,C,D) #no ":" at the end!
 	return ( (ccw(A, C, D) != ccw(B, C, D)) && ccw(A, B, C) != ccw(A, B, D))
end

function nearestN(r, nodeslist) #dear lord rewrite this so it doesn't need to pass in maxNodeId
    nearestDist = 9999;
    nearestNode = [];
    for n in nodeslist
        dist = distPt(r, n.state)
        if dist < nearestDist
            nearestDist = dist
            nearestNode = n
        end
    end
    return nearestNode
end

function distPt(pt1, pt2)
    x2,y2 = pt2.x, pt2.y
    x1,y1 = pt1.x, pt1.y
    dist = sqrt( (x1-x2)^2 + (y1-y2)^2 )
    return dist
end

function rrtPathPlanner(niterations)
    nIter = niterations
    maxDist = 3
    #room = Room(0,0,21,21);

    obs1 = rrt.Obstacle(rrt.Point(8,3),rrt.Point(10,18)) #Todo
    #obs2 = rrt.Obstacle(rrt.Point(10,12),rrt.Point(14,12)) #Todo
    obs1 = rrt.Obstacle(rrt.Point(8,5),rrt.Point(10,18)) #Todo

    rrtstart = rrt.Point(1,0)
    goal = rrt.Point(18,18)

    nodeslist = Vector{rrt.Node}()

    startNode = rrt.Node(0,0, rrtstart)
    push!(nodeslist, startNode)
    #@printf("string %s",nodeslist)

    maxNodeID = 1
    isPathFound = false

    for i in 1:nIter
        r = rrt.Point(rand(1:20),rand(1:20)) #new point
        nn = nearestN(r, nodeslist) # parent point

        #if !isCollidingNode(r, obs1)#check node XY first #Tofix: this function
        if !isCollidingNode(r, obs1) && !isCollidingEdge(r, nn.state, obs1) # check edge #rewrite so we can check multiple obstacles...
                nearestDist = distPt(r, nn.state)
                    if nearestDist > maxDist
                        x1,y1 = nn.state.x, nn.state.y
                        x2,y2 = r.x, r.y
                        tantheta = atan2((y2-y1) , (x2-x1))
                        newX = x1 + maxDist * cos(tantheta)
                        newY = y1 + maxDist * sin(tantheta)
                        newPt = rrt.Point(floor(newX), floor(newY)) #Todo: maybe I shouldn't floor?
                        node = rrt.Node(maxNodeID, nn.id, newPt) #throw out r, but try to steer toward it
                    else
                        node = rrt.Node(maxNodeID, nn.id, r)
                    end
                maxNodeID += 1
                push!(nodeslist, node)

                if inGoalRegion(node, goal)
                    #print("Goallllll! This is the winning node:\n")
                    #@show node  #winning node
                    #@printf("Goallll! Found after %d iterations", i)
                    isPathFound = true
                    break
                end
        else
            end
        end

    #@show nodeslist
    #@printf("\nPath found? %s Length of nodeslist: %d\n", isPathFound, length(nodeslist))
    cost = costWinningPath(nodeslist)
    return cost, isPathFound, nodeslist
end

function plotPath(isPathFound, nlist) #rewrite so don't need to pass in isPathFound, obs1, rrtstart, goal, room

    ### <COPIED FOR NOW #Todo fix hardcoding
    obs1 = rrt.Obstacle(rrt.Point(8,5),rrt.Point(10,18)) #Todo

    rrtstart = rrt.Point(1,0)
    goal = rrt.Point(18,18)
    ### / COPIED FOR NOW>

    foo = rand(1)
    h = plot()

    @printf("%s", "plotted\n")
    plot!(h, show=true, legend=false, size=(600,600),xaxis=((-5,25), 0:1:20 ), yaxis=((-5,25), 0:1:20), foreground_color_grid=:lightcyan)

    nIter = -1#fix hardcoding
    title!("RRT nIter = $(nIter), Path Found $(isPathFound)")

    # plot room
    dim = 21
    roomx = [0,0,dim,dim];
    roomy = [0,dim,dim,0];
    plot!(roomx, roomy, color=:black, linewidth=5)

    # plot start and end goals
    circle(1,0, 0.5, :red)
    circle(18,18, 1.4, :forestgreen)

    # plot obstacles
    circleObs(obs1)

    # plot all paths

    # plot winning path
    if isPathFound
        cost = plotWinningPath(nlist)
        @printf("\n!!!! the cost of the path was %d across %d nodes  !!!! \n", cost, length(nlist))
        #nEnd = nlist[end]
    end
    for n in nlist
        if n.id == 0
            continue
        else
            plotEdge(n,nlist)
        end
    end
    # display winning path cost
end

#module pHelp()

    ### Some helper functions
    function circle(x,y,r,c_color)
        th = 0:pi/50:2*pi;
        xunit = r * cos.(th) + x;
        yunit = r * sin.(th) + y;
        plot!(xunit, yunit,color=c_color,linewidth=3.0);
    end

    function circleObs(obstacle)
        x1,y1 = obstacle.SW.x, obstacle.SW.y
        x2,y2 = obstacle.NE.x, obstacle.NE.y
        r = 0.2
        obsColor = :blue

        circle(x1,y1,r,obsColor)
        circle(x1,y2,r,obsColor)
        circle(x2,y1,r,obsColor)
        circle(x2,y2,r,obsColor)

        plot!([x1,x1,x2,x2,x1], [y1,y2,y2,y1,y1], color=obsColor, linewidth=2)
    end

    function plotEdge(n, nlist)
        pt1 = n.state
        iPrev = n.iPrev
        nPrev = findNode(iPrev, nlist)
        pt2 = nPrev.state
        plot!( [pt1.x, pt2.x],[pt1.y, pt2.y], color=:orange, linewidth=3)
        # Todo: rewrite this to just plot them all at once.... [x1 x2 x3] [y1 y2 y3]
    end
    function costWinningPath(nlist)
        curNode = nlist[end]
        guhPath = [ rrt.Point(curNode.state.x, curNode.state.y)]
        while true
            iPrev = curNode.iPrev
            curNode = findNode(iPrev, nlist)
            x,y = curNode.state.x, curNode.state.y
            push!( guhPath, rrt.Point(x,y))
            if curNode.id == 0
                push!( guhPath, rrt.Point(x,y))
                break
            end
        end
        cost = 0
        for i in 2:length(guhPath)
            cost += distPt(guhPath[i],guhPath[i-1])
        end

        return cost
    end

    function plotWinningPath(nlist)
        curNode = nlist[end]
        xPath = [curNode.state.x]
        yPath = [curNode.state.y]
        guhPath = [ rrt.Point(curNode.state.x, curNode.state.y)]
        while true
            iPrev = curNode.iPrev
            curNode = findNode(iPrev, nlist)
            x,y = curNode.state.x, curNode.state.y
            push!(xPath, x)
            push!(yPath,y)
            push!( guhPath, rrt.Point(x,y))
            if curNode.id == 0
                push!(xPath, x)
                push!(yPath,y)
                break
            end
        end
        print("plotted winning path")
        cost = 0
        for i in 2:length(guhPath)
            cost += distPt(guhPath[i],guhPath[i-1])
        end

        plot!( xPath, yPath, color = :orchid, linewidth=3)
        #@show nlist
        return cost
    end

    function findNode(id, nodeslist)
        for n in nodeslist
            if n.id == id
                return n
            end
        end
    end
#end

cost, isPathFound, nlist = rrtPathPlanner(150) #maxIter
plotPath(isPathFound,nlist)


text:
#! /usr/bin/julia

# Rosetta Code, Polynomial long division

using Polynomials

p = Poly([-42,0,-12,1])
q = Poly([-3,1])

d, r = divrem(p,q)

println(p, " divided by ", q, " is ", d, " with remainder ", r, ".")

max_stars_repo_name: MichaeLeroy/rosetta-code

max_issues_repo_name: MichaeLeroy/rosetta-code

max_forks_repo_name: MichaeLeroy/rosetta-code


text:
"""
MomentUncertainty <: CenteredAmbiguitySet

```math
\\left\\{ r  \\; \\middle| \\begin{array}{ll}
s.t.  \\quad (\\mathbb{E} [r] - \\hat{r}) ' \\Sigma^{-1} (\\mathbb{E} [r] - \\hat{r}) \\leq \\gamma_1 \\\\
\\quad \\quad \\mathbb{E} [ (r - \\hat{r}) ' (r - \\hat{r}) ] \\leq \\gamma_2 \\Sigma \\\\
\\end{array}
\\right\\} \\\\
```

Atributes:
- `d::Sampleable{Multivariate, Continous}`: The parent distribution with an uncertain mean
- `γ1::Float64`: Uniform uncertainty around the mean (has to be greater than 0). (default: std(dist) / 5)
- `γ2::Float64`: Uncertainty around the covariance (has to be greater than 1). (default: 3.0)

References:
- Delage paper on moment uncertainty (implemented): https://www.researchgate.net/publication/220244490_Distributionally_Robust_Optimization_Under_Moment_Uncertainty_with_Application_to_Data-Driven_Problems

"""
struct MomentUncertainty{T<:Real, D<:ContinuousMultivariateSampleable} <: CenteredAmbiguitySet{T,D}
    d::D
    γ1::T
    γ2::T

    # Inner constructor for validating arguments
    function MomentUncertainty{T, D}(
        d::D, γ1::T, γ2::T
    ) where {T<:Real, D<:ContinuousMultivariateSampleable}
        γ1 >= 0 || throw(ArgumentError("γ1 must be >= 0"))
        γ2 >= 1 || throw(ArgumentError("γ2 must be >= 1"))
        return new{T, D}(d, γ1, γ2)
    end
end

# Default outer constructor
function MomentUncertainty(
    d::D, γ1::T, γ2::T
) where {T<:Real, D<:ContinuousMultivariateSampleable}
    MomentUncertainty{T, D}(d, γ1, γ2)
end

distribution(s::MomentUncertainty) = s.d

"""
    calculate_measure!(measure::ExpectedReturn{MomentUncertainty,WorstCase}, w)

Returns worst case utility return (WCR) under distribution uncertainty defined by MomentUncertainty ambiguity set ([`MomentUncertainty`](@ref)).

Arguments:
 - `w`: portfolio optimization investment variable ("weights").
 - `s::MomentUncertainty`: Struct containing atributes of MomentUncertainty ambiguity set.
"""
function calculate_measure!(m::ExpectedUtility{U,S,R}, w) where {U<:PieceWiseUtility,S<:MomentUncertainty,R<:WorstCase}
    model = owner_model(w)
    s = ambiguityset(m)
    utility_function = utility(m)

    # parameters
    means = Vector(Distributions.mean(s.d))
    Σ = Matrix(Distributions.cov(s.d))
    n = length(s)
    a = coefficients(utility_function)
    b = intercepts(utility_function)
    γ1 = s.γ1
    γ2 = s.γ2
    K = length(a)

    # dual variables
    @variable(model, P[i=1:n, j=1:n])
    @variable(model, p[i=1:n])
    @variable(model, s)
    @variable(model, Q[i=1:n, j=1:n])
    @variable(model, q[i=1:n])
    @variable(model, r)

    # constraints: from duality
    @constraint(model, p .== -q / 2 - Q * means)
    @constraint(model, [[P p]; [p' s]] in PSDCone())
    for k in 1:K
        @constraint(
            model, [[Q (q / 2 + a[k] * w / 2)]; [(q' / 2 + a[k] * w' / 2) (r + b[k])]] in PSDCone()
        )
    end

    return -(γ2 * dot(Σ, Q) - first(means'Q * means) + r + dot(Σ, P) - 2 * dot(means, p) + γ1 * s)
end

max_stars_repo_name: andrewrosemberg/PortfolioOpt

max_issues_repo_name: andrewrosemberg/PortfolioOpt

max_forks_repo_name: andrewrosemberg/PortfolioOpt


text:
@testset "$(rpad("Tableau coefficients",80))" begin

    using LinearAlgebra: normalize
    using GeometricBase
    using GeometricIntegrators.Utils

    import GeometricIntegrators.Tableaus: get_lobatto_nodes, get_lobatto_weights,
                                          get_lobatto_glrk_coefficients
    import QuadratureRules: LobattoLegendreQuadrature, nodes, weights

    function _get_lobatto_projective_stage(s, T=Float64)
        if s == 1
            a = reshape(T[0,1], (2,1))
        elseif s == 2
            a = T[ 0     0
                1//2  1//2 ]
        elseif s == 3
            a = T[ 0      0      0
                5//18  8//18  5//18 ]
        elseif s == 4
            a = T[ 0            0            0            0
                1//4-√30/72  1//4+√30/72  1//4+√30/72  1//4-√30/72 ]
        elseif s == 5
            a = T[ 0            0            0            0            0
                (-5797*√70 - 19635*√6 + 3297*√105 + 34069)/(900*(-31*√70 - 105*√6 + 12*√105 + 124))  3*(-155*√70 - 525*√6 + 21*√105 + 217)/(100*(-31*√70 - 105*√6 + 12*√105 + 124))  64/225  3*(-155*√70 - 525*√6 + 21*√105 + 217)/(100*(-31*√70 - 105*√6 + 12*√105 + 124))  (-5797*√70 - 19635*√6 + 3297*√105 + 34069)/(900*(-31*√70 - 105*√6 + 12*√105 + 124)) ]
        else
            @error "Lobatto projective coefficients for " * string(s) * " stages not implemented."
        end

        return a
    end

    @test get_lobatto_glrk_coefficients(1,2).a ≈ _get_lobatto_projective_stage(1)
    @test get_lobatto_glrk_coefficients(2,2).a ≈ _get_lobatto_projective_stage(2)
    @test get_lobatto_glrk_coefficients(3,2).a ≈ _get_lobatto_projective_stage(3)
    @test get_lobatto_glrk_coefficients(4,2).a ≈ _get_lobatto_projective_stage(4)
    @test get_lobatto_glrk_coefficients(5,2).a ≈ _get_lobatto_projective_stage(5)

    function _get_lobatto_interstage_coefficients(s, σ=s+1, T=Float64)
        if s == 1 && σ == 2
            a = reshape([
                    [0]
                    [1]
                ], σ, s)
        elseif s == 2 && σ == 3
            a = [
                    [0         0       ]
                    [1/4+√3/8  1/4-√3/8]
                    [1/2       1/2     ]
                ]
        elseif s == 3 && σ == 4
            a = [
                    [0                     0              0                 ]
                    [5/36-√5/180+√15/30    2/9-4*√5/45    5/36-√15/30-√5/180]
                    [5/36+√5/180+√15/30    2/9+4*√5/45    5/36+√5/180-√15/30]
                    [5/18                  4/9            5/18              ]
                ]
        elseif s == 4 && σ == 5
            a = [
                    [ 0  0  0  0 ]
                    [ 0.1591671294900345358547852036503968244037099929876660089701433805420016628303049  0.01565551541810189985834126767767053597372403599533987485875443079743596295789264  -0.002649931830622319490548503812025451589560291885905677307276227674989121280824921  0.0005004515684973118782758043605289134275222217051875147207182646279180365249293946 ]
                    [ 0.176105937938662021225385981377510389121610166295011341406551070207182545793024   0.3049159394838621526624258370570061164298248408582482915516281042450345033436905    0.02115663794741091865104218833199417995250081122480493820210723599558956292335403  -0.002178515369935092538854006766510685503935818378064571160286410447806612060067542  ]
                    [ 0.1734269710002296168082561702504707901901521262117592555255463951314578972080256  0.3287225092618953908040165292010257479718859439689589070610115679156131875478712    0.3104170620131711714551267577113297604086016160877133548949809094431881033091524    0.01476029307869239283174677096060287921396435492928076127612127921737427090265032   ]
                    [ 0.1739274225687269286865319746109997036176743479169467702462646597593759337329541  0.3260725774312730713134680253890002963823256520830532297537353402406240662670459    0.3260725774312730713134680253890002963823256520830532297537353402406240662670459    0.1739274225687269286865319746109997036176743479169467702462646597593759337329541    ]
                ]
        else
            @error("Number of stages s=$(s) and σ=$(σ) is not supported.")
        end

        CoefficientsIRK{T}(:LobattoIIIIS, s^2, s, σ, a, get_lobatto_weights(σ), get_lobatto_nodes(σ))
    end

    @test get_lobatto_glrk_coefficients(1) ≈ _get_lobatto_interstage_coefficients(1)
    @test get_lobatto_glrk_coefficients(2) ≈ _get_lobatto_interstage_coefficients(2)
    @test get_lobatto_glrk_coefficients(3) ≈ _get_lobatto_interstage_coefficients(3)
    @test get_lobatto_glrk_coefficients(4) ≈ _get_lobatto_interstage_coefficients(4)

    # test PGLRK coefficients

    @test typeof(CoefficientsPGLRK(2)) <: CoefficientsPGLRK

end


text:
"""
	normalizemean!(A)

In place removal of the variable mean from the P×N matrix `A` where P is the number of variables and N is the number of samples.
"""
normalizemean!(A::AbstractMatrix) = A .-= mean(A,dims=2)

"""
	normalizemean(A)

Remove the variable mean from the P×N matrix `A` where P is the number of variables and N is the number of samples.
"""
normalizemean(A::AbstractMatrix) = A .- mean(A,dims=2)

"""
	normalizemeanstd!(A)

In place normalization of variables to be mean zero and standard deviation one of the P×N matrix `A` where P is the number of variables and N is the number of samples.
"""
function normalizemeanstd!(A::AbstractMatrix)
	normalizemean!(A::AbstractMatrix)
	A ./= max.(1e-12, std(A,dims=2)) # max to avoid div by zero
end
"""
	normalizemeanstd!(A)

Normalize variables to be mean zero and standard deviation one in the P×N matrix `A` where P is the number of variables and N is the number of samples.
"""
function normalizemeanstd(A::AbstractMatrix)
	A = normalizemean(A)
	A ./= max.(1e-12, std(A,dims=2)) # max to avoid div by zero
end


text:
using TableView
using Test, WebIO, Tables, JSON

@testset "installation" begin
    version = readchomp(joinpath(@__DIR__, "..", "ag-grid.version"))
    @test isfile(joinpath(@__DIR__, "..", "deps", "ag-grid-$(version)", "ag-grid.js"))
end
@testset "named tuple table" begin
    nttable = [
        (a = 2.0, b = 3),
        (a = 3.0, b = 12)
    ]
    @test showtable(nttable) isa WebIO.Scope
end
@testset "named tuple table with missing and nothing" begin
    nttable = [
        (a = 2.0, b = 3, c = missing),
        (a = 3.0, b = 12, c = nothing)
    ]
    @test showtable(nttable) isa WebIO.Scope
end
@testset "inf and nan serializing" begin
    rows = Tables.table([NaN Inf -Inf 0])
    names = [:a, :b, :c, :d, :e, :f]
    types = [Float64 for _ in 1:6]
    json = TableView.table2json(rows, names, types)
    firstrow = JSON.parse(json)[1]
    @test firstrow["a"] == "NaN"
    @test firstrow["b"] == "Inf"
    @test firstrow["c"] == "-Inf"
    @test firstrow["d"] == 0
end
@testset "normal array" begin
    array = rand(10, 10)
    @test showtable(array) isa WebIO.Scope
end


text:
module FrostHelperRainbowTilesetController

using ..Ahorn, Maple

@mapdef Entity "FrostHelper/RainbowTilesetController" RainbowTilesetController(x::Integer, y::Integer, tilesets::String="3", bg::Bool=false)

const placements = Ahorn.PlacementDict(
	"Rainbow Tileset Controller (Frost Helper)" => Ahorn.EntityPlacement(
		RainbowTilesetController
	)
)

const sprite = "editor/FrostHelper/RainbowTilesetController"

function Ahorn.selection(entity::RainbowTilesetController)
	x, y = Ahorn.position(entity)

	return Ahorn.getSpriteRectangle(sprite, x, y)
end

function Ahorn.renderAbs(ctx::Ahorn.Cairo.CairoContext, entity::RainbowTilesetController, room::Maple.Room)
	x, y = Ahorn.position(entity)

	Ahorn.drawSprite(ctx, sprite, x, y)
end

end

max_stars_repo_name: JaThePlayer/FrostHelper

max_issues_repo_name: JaThePlayer/FrostHelper

max_forks_repo_name: JaThePlayer/FrostHelper


text:
using Resonance

metadata = CSV.read("data/wrangled.csv", DataFrame)
@rsubset!(metadata, !ismissing(:subject),
                    !ismissing(:sample),
                    )
@rsubset!(metadata, !startswith(:sample, "C"),
                    !startswith(:sample, "z"),
                    )

volumes = metadata[:, r"^(left|right|sample)"i][:, 2:end]
dropmissing!(volumes)

spec = metaphlan_profiles(filter(f-> contains(f, "profile"), readdir("/grace/echo/analysis/biobakery3/links/metaphlan", join=true)), :species)

sns = map(samplenames(spec)) do s
    replace(s, r"_S\d+_profile"=>"")
end

idx = unique([findfirst(==(x), sns) for x in sns])

spec = CommunityProfile(abundances(spec)[:, idx], features(spec), map(samplenames(spec)[idx]) do s
    MicrobiomeSample(replace(s, r"_S\d+_profile"=>""))
end)

overlap = sort(samplenames(spec) ∩ volumes.sample)

brainspec = spec[:, overlap]

volumes = @chain volumes begin
    @rsubset(:sample in overlap)
    @orderby(:sample)
end

permutedims(volumes, "sample")

CSV.write("data/halla_volumes.tsv", permutedims(volumes, :sample), delim='\t')
CSV.write("data/halla_species.tsv", brainspec, delim='\t')

##

run(`halla -x data/halla_species.tsv -y data/halla_volumes.tsv --alla -o data/halla_brain/`)

##

metabs = CSV.read("data/metabolites.csv", DataFrame)

##

sample_meta = airtable_metadata()
@rsubset!(sample_meta, in(:sid_old, names(metabs)))

rename!(metabs, Dict(row.sid_old => row.sample for row in eachrow(sample_meta)))
rename!(sample_meta, "sample"=>"etoh_sample")

sample_meta = leftjoin(sample_meta[:, [:subject, :timepoint, :etoh_sample]], metadata, on=[:subject,:timepoint])

@rsubset!(sample_meta, !ismissing(:sample) && in(:sample, samplenames(spec)))
metabspec = let cmp = spec[:, string.(sample_meta.sample)]
    CommunityProfile(abundances(cmp), features(cmp), MicrobiomeSample.(string.(sample_meta.etoh_sample)))
end

CSV.write("data/halla_metabolites.tsv", select(metabs, ["Metabolite", sample_meta.etoh_sample...]), delim='\t')
CSV.write("data/halla_metab_species.tsv", metabspec, delim='\t')

##

run(`halla -x data/halla_metab_species.tsv -y data/halla_metabolites.tsv --alla -o data/halla_metab/`)

##

using CairoMakie

fig = Figure(resolution=(1200,1200))

ax, sc = scatter(fig[1,1], volumes."Right-Cerebral-Cortex", volumes."Left-Cerebral-Cortex")

fig

CSV.write("data/volumes_desc.csv", describe(volumes))

max_stars_repo_name: Klepac-Ceraj-Lab/ResonanceAnalysis

max_issues_repo_name: Klepac-Ceraj-Lab/ResonanceAnalysis

max_forks_repo_name: Klepac-Ceraj-Lab/ResonanceAnalysis


text:
"""
    Corpus{T<:NamedTuple,D}

A corpus is a collection of [`Document`](@ref)s, along with some metadata. It has two fields,

* `documents::Vector{Document{D}}`
* `metadata::T`

Note each `Document` in a `Corpus` must have metadata of the same type.
"""
struct Corpus{T<:NamedTuple,D}
    documents::Vector{Document{D}}
    metadata::T

    function Corpus(documents::Vector{Document{D}}, metadata::T) where {T,D}
        disjoint_keys_check(T, D)
        check_keys(T)

        return new{T,D}(documents, metadata)
    end
end

function Base.match(query::AbstractQuery, p::Corpus)
    for r in p.documents
        m = match(query, r)
        if m !== nothing
            return m
        end
    end
    return nothing
end

function match_all(query::AbstractQuery, p::Corpus)
    return reduce(vcat, (match_all(query, document) for document in p.documents))
end

"""
    Base.match(query::AbstractQuery, corpus::Corpus)

Looks for a match for `query` in any [`Document`](@ref) in `corpus`.
Returns either `nothing` if no match is found in any `Document`,
or a [`QueryMatch`](@ref) object.
"""
Base.match(query::AbstractQuery, corpus::Corpus)

"""
    match_all(query::AbstractQuery, corpus::Corpus)

Looks for all matches for `query` from all documents in `corpus`. Returns a
`Vector` of `QueryMatch` objects corresponding to all of the matches found,
across all doucments.
"""
match_all(query::AbstractQuery, corpus::Corpus)

max_stars_repo_name: beacon-biosignals/KeywordSearch

max_issues_repo_name: beacon-biosignals/KeywordSearch

max_forks_repo_name: beacon-biosignals/KeywordSearch


text:
# ---
# title: Discrete scales
# cover: assets/discrete_scales.png
# description: Sorting and renaming categorical variables
# author: "[Pietro Vertechi](https://github.com/piever)"
# ---

using AlgebraOfGraphics, CairoMakie
set_aog_theme!() #src

# By default categorical ticks, as well as names from legend entries, are taken from the
# value of the variable converted to a string. Scales can be equipped with labels to
# overwrite that

df = (x=rand(["a", "b", "c"], 100), y=rand(100))
plt = data(df) * mapping(:x, :y) * visual(BoxPlot)
draw(plt)

#

df = (x=rand(["a", "b", "c"], 100), y=rand(100))
plt = data(df) *
    mapping(
        :x => renamer("a" => "label1", "b" => "label2", "c" => "label3"),
        :y
    ) * visual(BoxPlot)
draw(plt)

# The order can also be changed by tweaking the scale

plt = data(df) *
    mapping(
        :x => renamer("b" => "label b", "a" => "label a", "c" => "label c"),
        :y
    ) * visual(BoxPlot)
fg = draw(plt)

# save cover image #src
mkpath("assets") #src
save("assets/discrete_scales.png", fg) #src


text:
InputFolder = './Images/VesselsSpots/';
OutputFolder = './Results/Images/VesselsSpots/Vessels/';
Fill = -1;
ExportDist = 1;

@iA = '*.tif';

@fxm_mFilterObjSize [iA] > [M];
params.MinArea = 2;
params.MaxArea = Inf;
/endf

/show iA > M;
/keep M > tif;

max_stars_repo_name: SebastienTs/LOBSTER

max_issues_repo_name: SebastienTs/LOBSTER

max_forks_repo_name: SebastienTs/LOBSTER