text:
"""
    DeforModelRedModule

Module for  model reduction in models of deformation.
"""
module DeforModelRedModule

using FinEtools.FTypesModule: FInt, FFlt, FCplxFlt, FFltVec, FIntVec, FFltMat, FIntMat, FMat, FVec, FDataDict

abstract type DeforModelRed end;
struct DeforModelRed1D <: DeforModelRed
end
struct DeforModelRed2DStrain <: DeforModelRed
end
struct DeforModelRed2DStress <: DeforModelRed
end
struct DeforModelRed2DAxisymm <: DeforModelRed
end
struct DeforModelRed3D <: DeforModelRed
end

"""
    nstressstrain(::Type{DeforModelRed})::FInt

How many stress components involved in the balance equation (and how many
strain components) are there for each model reduction type?
"""
nstressstrain(::Type{DeforModelRed1D})::FInt = 1;
nstressstrain(::Type{DeforModelRed2DStrain})::FInt = 3;
nstressstrain(::Type{DeforModelRed2DStress})::FInt = 3;
nstressstrain(::Type{DeforModelRed2DAxisymm})::FInt = 4;
nstressstrain(::Type{DeforModelRed3D})::FInt = 6;

"""
    nthermstrain(::Type{DeforModelRed1D})::FInt

How many thermal strain components are there for each model reduction type?
"""
nthermstrain(::Type{DeforModelRed1D})::FInt = 1;
nthermstrain(::Type{DeforModelRed2DStrain})::FInt = 4;
nthermstrain(::Type{DeforModelRed2DStress})::FInt = 3;
nthermstrain(::Type{DeforModelRed2DAxisymm})::FInt = 4;
nthermstrain(::Type{DeforModelRed3D})::FInt = 6;

"""
    stresscomponentmap(::Type{DeforModelRed1D})

Construct a dictionary to map from stress-component symbols to indexes.
"""
function stresscomponentmap(::Type{DeforModelRed1D})
    a=Dict{Symbol,FInt}();
    a[:x]=1;
    a[:xx]=1;
    return a
end

function stresscomponentmap(::Type{DeforModelRed2DStrain})
    a=Dict{Symbol,FInt}();
    a[:x]=1;  a[:y]=2;  a[:xy]=3;  a[:z]=4;
    a[:xx]=1;  a[:yy]=2; a[:zz]=4;
    return a
end

function stresscomponentmap(::Type{DeforModelRed2DStress})
    a=Dict{Symbol,FInt}();
    a[:x]=1;  a[:y]=2;  a[:xy]=3;
    a[:xx]=1;  a[:yy]=2;
    return a
end

function stresscomponentmap(::Type{DeforModelRed2DAxisymm})
    a=Dict{Symbol,FInt}();
    a[:x]=1;  a[:y]=2;  a[:z]=3;  a[:xy]=4;
    a[:xx]=1;  a[:yy]=2;  a[:zz]=3;
    return a
end

function stresscomponentmap(::Type{DeforModelRed3D})
    a=Dict{Symbol,FInt}();
    a[:x]=1;  a[:y]=2;  a[:z]=3;  a[:xy]=4; a[:xz]=5; a[:yz]=6;
    a[:xx]=1;  a[:yy]=2;  a[:zz]=3;
    return a
end

"""
    Blmat!(MR::Type{DeforModelRed1D}, B::FFltMat,
      N::FFltMat, gradN::FFltMat, c::FFltMat, Rm::FFltMat)

Compute the strain-displacement matrix for a one-manifold element.

Compute the linear, displacement independent, strain-displacement matrix for a
one-manifold element.   *The input displacements are in the global Cartesian
coordinate system, the output strains are in the material coordinate system.*

`N` =matrix of basis function values
`gradN` =matrix of basis function gradients with respect to the
  Cartesian coordinates in the directions of the material orientation
`c` =array of spatial coordinates of the evaluation point
  in the global Cartesian coordinates.
`Rm` =orthogonal matrix with the unit basis vectors of the local
  material orientation coordinate system as columns.
  `size(Rm)= [ndim,mdim]`, where `ndim` = number of spatial dimensions
  of the embedding space (here `ndim >= 3`), aand `mdim` = number of manifold
  dimensions (here `mdim == 1`).

Output:
`B` = strain-displacement matrix, where  `size(B) = [nstressstrain,nnodes*dim]`;
  here nstressstrain= number of strains, dim = Number of spatial
  dimensions of the embedding space, and nnodes = number of
  finite element nodes on the element. The strain components are ordered as
  shown  in `stresscomponentmap`. The matrix is passed in as a buffer,
  set to zero,  and filled in  with the nonzero components.  It is also returned
for convenience.
"""
function Blmat!(MR::Type{DeforModelRed1D}, B::FFltMat,
  N::FFltMat, gradN::FFltMat, c::FFltMat, Rm::FFltMat)
  nnodes = size(gradN,1);
  @assert (dim = size(c,2)) >= 1;
  @assert size(B) == (nstressstrain(MR), dim*nnodes)
  fill!(B,0.0);
  @inbounds for i= 1:nnodes
    @inbounds for j=1:dim
      k=dim*(i-1);
      B[1,k+j]=  gradN[i,1] *Rm[j,1];
    end
  end
  return B::FFltMat;
end

"""
    Blmat!(MR::Type{DeforModelRed2DStrain}, B::FFltMat,
      N::FFltMat, gradN::FFltMat, c::FFltMat, Rm::FFltMat)

Compute the strain-displacement matrix for a two-manifold element for
plane-strain models.

Compute the linear, displacement independent, strain-displacement matrix for a
two-manifold element in plane-strain setting.   *The input displacements
are in the global Cartesian coordinate system, the output strains are in the
material coordinate system.*

`N` =matrix of basis function values
`gradN` =matrix of basis function gradients with respect to the
  Cartesian coordinates in the directions of the material orientation
`c` =array of spatial coordinates of the evaluation point
  in the global Cartesian coordinates.
`Rm` =orthogonal matrix with the unit basis vectors of the local
  material orientation coordinate system as columns.
  `size(Rm)= [ndim,mdim]`, where `ndim` = number of spatial dimensions
  of the embedding space (here `ndim == 2`), aand `mdim` = number of manifold
  dimensions (here `mdim == 2`).

Output:
`B` = strain-displacement matrix, where  `size(B) = [nstressstrain,nnodes*dim]`;
  here nstressstrain= number of strains, dim = Number of spatial
  dimensions of the embedding space, and nnodes = number of
  finite element nodes on the element. The strain components are ordered as
  shown  in `stresscomponentmap`. The matrix is passed in as a buffer,
  set to zero,  and filled in  with the nonzero components.  It is also returned
for convenience.
"""
function Blmat!(MR::Type{DeforModelRed2DStrain}, B::FFltMat,
  N::FFltMat, gradN::FFltMat, c::FFltMat, Rm::FFltMat)
  nnodes = size(gradN,1);
  @assert (dim = size(c,2)) == 2;
  @assert size(B) == (nstressstrain(MR), dim*nnodes)
  fill!(B,0.0);
  @inbounds for i= 1:nnodes
    k=dim*(i-1);
    @inbounds for j=1:dim
      B[1,k+j]=gradN[i,1]*Rm[j,1]
      B[2,k+j]=gradN[i,2]*Rm[j,2]
      B[3,k+j]=gradN[i,2]*Rm[j,1]+gradN[i,1]*Rm[j,2]
    end
  end
  return B::FFltMat;
end

"""
    Blmat!(MR::Type{DeforModelRed2DStress}, B::FFltMat,
      N::FFltMat, gradN::FFltMat, c::FFltMat, Rm::FFltMat)

Compute the strain-displacement matrix for a two-manifold element for
plane-stress models.

Compute the linear, displacement independent, strain-displacement matrix for a
two-manifold element in plane-stress setting.   *The input displacements
are in the global Cartesian coordinate system, the output strains are in the
material coordinate system.*

`N` =matrix of basis function values
`gradN` =matrix of basis function gradients with respect to the
  Cartesian coordinates in the directions of the material orientation
`c` =array of spatial coordinates of the evaluation point
  in the global Cartesian coordinates.
`Rm` =orthogonal matrix with the unit basis vectors of the local
  material orientation coordinate system as columns.
  `size(Rm)= [ndim,mdim]`, where `ndim` = number of spatial dimensions
  of the embedding space (here `ndim >= 3`), aand `mdim` = number of manifold
  dimensions (here `mdim=2`).

Output:
`B` = strain-displacement matrix, where  `size(B) = [nstressstrain,nnodes*dim]`;
  here nstressstrain= number of strains, dim = Number of spatial
  dimensions of the embedding space, and nnodes = number of
  finite element nodes on the element. The strain components are ordered as
  shown  in `stresscomponentmap`. The matrix is passed in as a buffer,
  set to zero,  and filled in  with the nonzero components.  It is also returned
for convenience.
"""
function Blmat!(MR::Type{DeforModelRed2DStress}, B::FFltMat,
  N::FFltMat, gradN::FFltMat, c::FFltMat, Rm::FFltMat)
  nnodes = size(gradN,1);
  @assert (dim = size(c,2)) >= 2;
  @assert size(B) == (nstressstrain(MR), dim*nnodes)
  fill!(B,0.0);
  @inbounds for i= 1:nnodes
    k=dim*(i-1);
    @inbounds for j=1:dim
      B[1,k+j]=gradN[i,1]*Rm[j,1]
      B[2,k+j]=gradN[i,2]*Rm[j,2]
      B[3,k+j]=gradN[i,2]*Rm[j,1]+gradN[i,1]*Rm[j,2]
    end
  end
  return B::FFltMat;
end

"""
    Blmat!(MR::Type{DeforModelRed2DAxisymm}, B::FFltMat,
      N::FFltMat, gradN::FFltMat, c::FFltMat, Rm::FFltMat)

Compute the strain-displacement matrix for a two-manifold element for
axially symmetric models.

Compute the linear, displacement independent, strain-displacement matrix for a
two-manifold element in axially-symmetric setting.   *The input displacements
are in the global Cartesian coordinate system, the output strains are in the
material coordinate system.*

`N` =matrix of basis function values
`gradN` =matrix of basis function gradients with respect to the
  Cartesian coordinates in the directions of the material orientation
`c` =array of spatial coordinates of the evaluation point
  in the global Cartesian coordinates; `c[1]` is the distance from the axis of
  symmetry.
`Rm` =orthogonal matrix with the unit basis vectors of the local
  material orientation coordinate system as columns.
  `size(Rm)= [ndim,ndim]`, where `ndim` = number of spatial dimensions
  of the embedding space (here `ndim=3`).

Output:
`B` = strain-displacement matrix, where  `size(B) = [nstressstrain,nnodes*dim]`;
  here nstressstrain= number of strains, dim = Number of spatial
  dimensions of the embedding space, and nnodes = number of
  finite element nodes on the element. The strain components are ordered as
  shown  in `stresscomponentmap`. The matrix is passed in as a buffer,
  set to zero,  and filled in  with the nonzero components.  It is also returned
for convenience.
"""
function Blmat!(MR::Type{DeforModelRed2DAxisymm}, B::FFltMat,
  N::FFltMat, gradN::FFltMat, c::FFltMat, Rm::FFltMat)
  nnodes = size(gradN,1);
  @assert (dim = size(c,2)) == 2;
  @assert size(B) == (nstressstrain(MR), dim*nnodes)
  r = c[1]; # distance from the axis of symmetry
  # if practically zero, clamped it to machine epsilon because we need to divide with it
  if r <= eps(1.0)
    r = eps(1.0);
  end
  fill!(B,0.0);
  @inbounds for i= 1:nnodes
    k=dim*(i-1);
    @inbounds for j=1:dim
      B[1,k+j]=gradN[i,1]*Rm[j,1]
      B[2,k+j]=gradN[i,2]*Rm[j,2]
      B[3,k+j]=N[i]/r*Rm[j,1]
      B[4,k+j]=gradN[i,2]*Rm[j,1]+gradN[i,1]*Rm[j,2]
    end
  end
  return B::FFltMat;
end

"""
    Blmat!(MR::Type{DeforModelRed3D}, B::FFltMat,
      N::FFltMat, gradN::FFltMat, c::FFltMat, Rm::FFltMat)

Compute the strain-displacement matrix for a three-manifold element.

Compute the linear, displacement independent, strain-displacement matrix
for a three-manifold element.   *The input displacements are in the global
Cartesian coordinate system, the output strains are in the material coordinate
system.*

`N` =matrix of basis function values
`gradN` =matrix of basis function gradients with respect to the
  Cartesian coordinates in the directions of the material orientation
`c` =array of spatial coordinates of the evaluation point
  in the global Cartesian coordinates
`Rm` =orthogonal matrix with the unit basis vectors of the local
  material orientation coordinate system as columns.
  `size(Rm)= [ndim,ndim]`, where `ndim` = number of spatial dimensions
  of the embedding space (here `ndim=3`).

Output:
`B` = strain-displacement matrix, where  `size(B) = [nstressstrain,nnodes*dim]`;
  here nstressstrain= number of strains, dim = Number of spatial
  dimensions of the embedding space, and nnodes = number of
  finite element nodes on the element. The strain components are ordered as
  shown  in `stresscomponentmap`. The matrix is passed in as a buffer,
  set to zero,  and filled in  with the nonzero components.  It is also returned
for convenience.
"""
function Blmat!(MR::Type{DeforModelRed3D}, B::FFltMat,
  N::FFltMat, gradN::FFltMat, c::FFltMat, Rm::FFltMat)
  nnodes = size(gradN,1);
  @assert (dim = size(c,2)) == 3;
  @assert size(B) == (nstressstrain(MR), dim*nnodes)
  fill!(B,0.0);
  @inbounds for i= 1:nnodes
    k=dim*(i-1);
    @inbounds for j=1:dim
      kj=k+j
      B[1,kj]=gradN[i,1]*Rm[j,1]
      B[2,kj]=gradN[i,2]*Rm[j,2]
      B[3,kj]=gradN[i,3]*Rm[j,3]
      B[4,kj]=gradN[i,2]*Rm[j,1]+gradN[i,1]*Rm[j,2]
      B[5,kj]=gradN[i,3]*Rm[j,1]+gradN[i,1]*Rm[j,3]
      B[6,kj]=gradN[i,3]*Rm[j,2]+gradN[i,2]*Rm[j,3]
    end
  end
  return B::FFltMat;
end

end


text:
const _RAND_STR_ALPB = ['a':'z'; 'A':'Z'; '0':'9']

rand_str(n = 10) = join(rand(_RAND_STR_ALPB, n))

ext: jl

max_stars_repo_name: josePereiro/GitWorkers

max_issues_repo_name: josePereiro/GitWorkers

max_forks_repo_name: josePereiro/GitWorkers


text:
# ---
# title: 27. Remove Element
# id: problem27
# author: zhwang
# date: 2022-03-16
# difficulty: Easy
# categories: Array, Two Pointers
# link: <https://leetcode.com/problems/remove-element/description/>
# hidden: true
# ---
#
# Given an array _nums_ and a value _val_ , remove all instances of that value
# [**in-place**](https://en.wikipedia.org/wiki/In-place_algorithm) and return
# the new length.
#
# Do not allocate extra space for another array, you must do this by **modifying
# the input array[in-place](https://en.wikipedia.org/wiki/In-place_algorithm)**
# with O(1) extra memory.
#
# The order of elements can be changed. It doesn't matter what you leave beyond
# the new length.
#
# **Clarification:**
#
# Confused why the returned value is an integer but your answer is an array?
#
# Note that the input array is passed in by **reference** , which means a
# modification to the input array will be known to the caller as well.
#
# Internally you can think of this:
#
#
#
#     // **nums** is passed in by reference. (i.e., without making a copy)
#     int len = removeElement(nums, val);
#
#     // any modification to **nums** in your function would be known by the caller.
#     // using the length returned by your function, it prints the first **len** elements.
#     for (int i = 0; i < len; i++) {
#         print(nums[i]);
#     }
#
#
#
# **Example 1:**
#
#
#
#     Input: nums = [3,2,2,3], val = 3
#     Output: 2, nums = [2,2]
#     Explanation: Your function should return length = **2** , with the first two elements of _nums_ being **2**.
#     It doesn't matter what you leave beyond the returned length. For example if you return 2 with nums = [2,2,3,3] or nums = [2,3,0,0], your answer will be accepted.
#
#
# **Example 2:**
#
#
#
#     Input: nums = [0,1,2,2,3,0,4,2], val = 2
#     Output: 5, nums = [0,1,4,0,3]
#     Explanation: Your function should return length = **5** , with the first five elements of _nums_ containing  **0** , **1** , **3** , **0** , and  **4**. Note that the order of those five elements can be arbitrary. It doesn't matter what values are set beyond the returned length.
#
#
#
#
# **Constraints:**
#
#   * `0 <= nums.length <= 100`
#   * `0 <= nums[i] <= 50`
#   * `0 <= val <= 100`
#
#
## @lc code=start
using LeetCode

remove_element!(nums::Vector{Int}, val::Int) = length(filter!(!=(val), nums))
## @lc code=end


text:
using EcoTrialStructure
using OffsetArrays
using Test
using Documenter

@testset "EcoTrialStructure.jl" begin
    offsetrange((indices, values)) = OffsetArrays.IdOffsetRange(;values, indices)

    @test isempty(detect_ambiguities(EcoTrialStructure))
    doctest(EcoTrialStructure)

    @testset "CellsTrial" begin
        function checkconsistent(ct, fs, cidx; errors::Bool=false)
            if errors
                @test !checkbounds(Bool, ct, fs, cidx)
                @test_throws BoundsError ct[fs,cidx]
            else
                @test  checkbounds(Bool, ct, fs, cidx)
                @test ct[fs,cidx][2] isa AbstractVecOrMat
            end
        end

        t = (100:100:500) * ms   # times at which frames were collected
        dFoF = [ 0.1 0.8;
                -0.1 0.7;
                0.2 0.6;
                0.1 0.5;
                0.0 0.4]
        ct = @inferred(CellsTrial(t, dFoF))
        ct2 = @inferred(CellsTrial([100,200,300,400,500] * ms, dFoF))
        @test_throws ArgumentError("frame times must be increasing") CellsTrial([100,300,200,400,500] * ms, dFoF)
        @test_throws DimensionMismatch CellsTrial((100:100:600) * ms, dFoF)
        @test_throws DimensionMismatch CellsTrial((100:100:400) * ms, dFoF)
        @test sprint(show, ct) == "2 cells with 5 timepoints"
        # equality and hashing
        @test ct == ct2
        d = Dict(ct => 1)
        d[ct2] = 2
        @test length(d) == 1 && d[ct] == 2
        # indexing
        @test ncells(ct) == 2
        @test ct[2:3,:] == (t[2:3], dFoF[2:3,:])
        @test ct[2:3,1] == (t[2:3], dFoF[2:3,1])
        @test ct[2:3,2] == (t[2:3], dFoF[2:3,2])
        @test ct[:,:]   == (t, dFoF)
        @test_throws BoundsError ct[2:3,3]
        @test_throws BoundsError ct[2:3,0]
        @test_throws BoundsError ct[4:6,:]
        @test FrameSeq(100ms, 2) == FrameSeq(100ms, 0:1)
        @test FrameSeq(:go, 2) == FrameSeq(:go, 0:1)
        fs = FrameSeq(:go, 2)
        @test length(fs) == 2
        fs = FrameSeq(100ms, 2)
        @test axes(ct[FrameSeq(100ms, 2),:][2], 1) == axes(fs)[1] == axes(fs, 1)
        @test ct[FrameSeq(100ms, 2),:] == ct[FrameSeq(120ms, 2),:] == (100ms..200ms, dFoF[offsetrange(0:1 => 1:2),:])
        @test ct[FrameSeq(200ms, -1:0),:] == (100ms..200ms, dFoF[offsetrange(-1:0 => 1:2),:])
        @test ct[FrameSeq(200ms, 2),:] == ct[FrameSeq(220ms, 2),:] == ct[FrameSeq(180ms, 2),:] == (200ms..300ms, dFoF[offsetrange(0:1 => 2:3),:])
        checkconsistent(ct, FrameSeq(100ms, 2), :)
        checkconsistent(ct, FrameSeq( 80ms, 2), :; errors=true)
        checkconsistent(ct, FrameSeq(100ms, -1:0), :; errors=true)
        checkconsistent(ct, FrameSeq(100ms, 3:4), :)
        checkconsistent(ct, FrameSeq(100ms, 4:5), :; errors=true)
        checkconsistent(ct, FrameSeq(470ms, 2), :; errors=true)
        checkconsistent(ct, FrameSeq(420ms, 2), :)
        @test ct[FrameSeq(490ms, 1),:] == (500ms..500ms, OffsetArray(dFoF[end:end,:], 0:0, Base.OneTo(size(dFoF,2))))
        @test_throws ArgumentError("indexing requires a concrete `FrameSeq`, use `fs(et::EventTiming)`") ct[FrameSeq(:go, 4), :]
    end

    @testset "OfferType & TrialResult" begin
        tt = OfferType(3, 1, false)
        @test OfferType(tt) === tt
        @test sprint(show, tt) == "OfferType(nA=3, nB=1, leftA=false)"
        @test eval(Meta.parse("OfferType(nA=3, nB=1, leftA=false)")) == tt

        tr = TrialResult(3, 1, false, true)
        @test TrialResult(tr) === tr
        @test OfferType(tr) == tt
        @test sprint(show, tr) == "TrialResult(nA=3, nB=1, leftA=false, choseA=true)"
        @test eval(Meta.parse("TrialResult(nA=3, nB=1, leftA=false, choseA=true)")) == tr
        @test tr == TrialResult(choseA=true, nB=1, nA=3, leftA=false)
        @test tr != TrialResult(3, 1, false, false)
        @test tr != TrialResult(3, 1, false, missing)
        @test TrialResult(3, 1, false, missing) == TrialResult(3, 1, false, missing)
        @test !isforced(tr)
        @test  madechoice(tr)
        @test !madechoice(TrialResult(3, 1, false, missing))
        ttf = TrialResult(3, 0, false, true)
        @test  isforced(ttf)
        @test !iswrong(ttf)
        ttf = TrialResult(3, 0, false, false)
        @test  isforced(ttf)
        @test  iswrong(ttf)
        @test !isforced(TrialResult(0, 0, false, true))

        # Ordering: useful for printing results along a 1d axis
        @test !(OfferType(5, 0, false) < OfferType(5, 0, false))
        @test !(OfferType(5, 0, false) > OfferType(5, 0, false))
        @test   OfferType(5, 0, false) < OfferType(4, 0, false)
        @test !(OfferType(5, 0, false) > OfferType(4, 0, false))
        @test   OfferType(5, 0, false) < OfferType(5, 2, false)
        @test   OfferType(5, 1, false) > OfferType(5, 0, false)

        @test   OfferType(5, 1, false) < OfferType(5, 2, false)
        @test !(OfferType(5, 1, false) > OfferType(5, 2, false))
        @test !(OfferType(5, 1, false) < OfferType(5, 1, false))
        @test !(OfferType(5, 1, false) > OfferType(5, 1, false))

        @test !(OfferType(0, 5, false) < OfferType(0, 5, false))
        @test !(OfferType(0, 5, false) > OfferType(0, 5, false))
        @test   OfferType(0, 5, false) > OfferType(0, 4, false)
        @test !(OfferType(0, 5, false) < OfferType(0, 4, false))
        @test   OfferType(0, 5, false) > OfferType(2, 5, false)
        @test   OfferType(1, 5, false) < OfferType(0, 5, false)

        @test !(OfferType(1, 5, false) < OfferType(1, 5, true))
        @test   OfferType(1, 5, true)  < OfferType(1, 5, false)
        @test   OfferType(2, 5, false) < OfferType(1, 5, true)
        @test   OfferType(2, 5, true)  < OfferType(1, 5, false)

        for (a, b) in ((true, true), (true, false), (false, true), (false, false), (true, missing), (missing, true))
            @test TrialResult(5, 1, false, a) < TrialResult(5, 2, false, b)    # OfferType trumps choice
        end
        tt = OfferType(5, 1, false)
        @test TrialResult(tt, true) < TrialResult(tt, missing) < TrialResult(tt, false)
        @test !(TrialResult(tt, true) < TrialResult(tt, true))
    end

    @testset "EventTiming" begin
        et = @inferred(EventTiming(0ms, 100ms, 400ms, 450ms, 837ms, 1.2s))
        @test sprint(show, et) == "EventTiming(trial_start=0.0f0 ms, offer_on=100.0f0 ms, offer_off=400.0f0 ms, go=450.0f0 ms, choice=837.0f0 ms, trial_end=1200.0f0 ms)"
        @test eval(Meta.parse("EventTiming(trial_start=0.0f0ms, offer_on=100.0f0ms, offer_off=400.0f0ms, go=450.0f0ms, choice=837.0f0ms, trial_end=1200.0f0ms)")) == et
        # Unitful doesn't yet support round-trip printing, see https://github.com/PainterQubits/Unitful.jl/pull/470
        @test_broken eval(Meta.parse("EventTiming(trial_start=0.0f0 ms, offer_on=100.0f0 ms, offer_off=400.0f0 ms, go=450.0f0 ms, choice=837.0f0 ms, trial_end=1200.0f0 ms)")) == et
        @test et.trial_end == 1200ms

        fs = FrameSeq(:go, 5)
        @test isdeferred(fs)
        fse = fs(et)
        @test fse.start == et.go
        @test fse(EventTiming(0s, 0s, 0s, 0s, 0s, 0s)) == fse   # if the timing is concrete, this doesn't change it
    end

    @testset "Matlab import" begin
        # Check the "reverse-encoding" of the TrialResult
        @test  EcoTrialStructure.mat_trialresult(1, 1, 1, 1).choseA
        @test !EcoTrialStructure.mat_trialresult(1, 1, 2, 1).choseA
        tr = EcoTrialStructure.mat_trialresult(0, 1, 0, 1)
        @test tr.choseA
        @test iswrong(tr)
        @test madechoice(tr)
        tr = EcoTrialStructure.mat_trialresult(1, 0, 0, 1)
        @test !tr.choseA
        @test iswrong(tr)
        @test madechoice(tr)
        @test_throws Exception EcoTrialStructure.mat_trialresult(1, 1, 0, 1)
        @test_throws Exception EcoTrialStructure.mat_trialresult(0, 0, 0, 1)
        tr = EcoTrialStructure.mat_trialresult(1, 1, 0, 0)
        @test !madechoice(tr)

        testfile = joinpath(@__DIR__, "data", "testfile.mat")  # snippet from data collected by Manning Zhang, Washington University in St. Louis
        cts, trs, ets = parsemat(testfile)
        tt3, tt201 = trs[3], trs[201]
        @test isforced(tt3)
        @test iswrong(tt3)
        @test tt3.nB == 2
        @test isforced(tt201)
        @test !iswrong(tt201)
        @test tt201.nA == 2
        @test sprint(show, cts[3])   == "63 cells with 69 timepoints"
        @test sprint(show, cts[201]) == "63 cells with 39 timepoints"
        @test length(ets) == 2
        et3 = ets[3]
        @test et3.trial_start == 0s
        @test et3.offer_on  == 2003ms
        @test et3.offer_off == et3.go == 4820ms
        @test et3.choice == 5111ms
        @test et3.trial_end == 11066ms

        cts, trs, ets, poscells = parsemat(positive_cells, testfile)
        @test sum(poscells) == 5 && length(poscells) == 63
    end
end


text:
using Colors
using FileIO

# All these palettes have 256 colors to simplify saving as 8-bit palettized PNG. (Maybe...)

"""
A colorblind friendly palette suited to use with log base-10 transformed data on the range [0.0, 1.0].
"""
const Log3BandColorblind = reduce(
    append!,
    (
        [colorant"rgb(0,0,0)"],
        range(colorant"rgb(0%,38%,100%)", stop = colorant"white", length = 120)[84:-1:1],
        range(colorant"rgb(75%,0%,0%)", stop = colorant"white", length = 120)[85:-1:1],
        range(colorant"rgb(100%,90%,0%)", stop = colorant"black", length = 85),
        [colorant"yellow"],
    ),
)

"""
Dave Bright's palette suited to use with log base-10 transformed data on the range [0.0, 1.0].
The 256ᵗʰ entry is yellow for an error condition.
"""
const Log3BandBright = reduce(
    append!,
    (
        [colorant"rgb(0,0,0)"],
        range(colorant"rgb(0,0,60)", stop = colorant"rgb(200,200,255)", length = 84), # 0–200/0–200/60–255
        range(colorant"rgb(0,60,0)", stop = colorant"rgb(200,255,200)", length = 85), # 0–200/60–255/0–200
        range(colorant"rgb(60,0,0)", stop = colorant"rgb(255,200,200)", length = 85), # 60–255/0–200/0–200
        [colorant"yellow"],
    ),
)

"""
An RGB gray-scale palette with entry 256 being yellow.
"""
const GrayScale = reduce(
    append!,
    (
        range(colorant"rgb(0,0,0)", stop = colorant"rgb(255,255,255)", length = 255),
        [colorant"yellow"],
    ),
)

"""
Tranforms numbers on the range [1.0e-3, 1.0] into a colorblind friendly palette using a log base-10 transform.
NaNs are plotted in yellow.
"""
Log3BandC(f::AbstractFloat)::Colorant = Log3BandColorblind[isnan(f) ? 256 :
                   min(
    max(trunc(Int, 255.0 * (log(10.0, max(f, 1.0e-3)) + 3.0) / 3.0), 0),
    254,
) + 1]

"""
Tranforms numbers on the range [1.0e-3, 1.0] into David Bright's Log3-band palette using a log base-10 transform.
NaNs are plotted in yellow.
"""
Log3Band(f::AbstractFloat)::Colorant = Log3BandBright[isnan(f) ? 256 :
               min(
    max(trunc(Int, 255.0 * (log(10.0, max(f, 1.0e-3)) + 3.0) / 3.0), 0),
    254,
) + 1]

"""
Tranforms numbers on the range [1.0e-3, 1.0] onto a Log base-10 gray scale palette. NaNs are plotted in yellow.
"""
LogScale(f::AbstractFloat)::Colorant = GrayScale[isnan(f) ? 256 :
          min(max(trunc(Int, 255.0 * (log(10.0, max(f, 1.0e-3)) + 3.0) / 3.0), 0), 254) + 1]

"""
Tranforms numbers on the range [0.0, 1.0] onto a linear gray scale palette. NaNs are plotted in yellow.
"""
LinearScale(f::AbstractFloat)::Colorant =
    GrayScale[isnan(f) ? 256 : min(254, max(0, trunc(Int, 254.0 * f))) + 1]

"""
Loads a legend image from the package source directory. Some legends include "LinearScale.png", "LogScale.png",
"Log3BandBright.png", "Log3BancColorblind.png".
"""
loadlegend(fn::String) = FileIO.load(joinpath(dirname(pathof(@__MODULE__)), fn))


text:
using BenchmarkTools

include("src/gauss_jordan_elim.jl");

# regular matrix
A = [1 2 3
     4 5 6];
gauss_jordan(A)

# a singular matrix
A = [1 2 3
     4 8 12];
gauss_jordan(A)

# less performatic because convert to float to avoid InexactError: Int64()
A = [1 4 8  1  4
     4 5 6  8  11
     1 3 2  4  8
     4 5 67 23 0];
@btime gauss_jordan(A)

# more performatic because it's already Float64
A = Float64[1 4 8  1  4
            4 5 6  8  11
            1 3 2  4  8
            4 5 67 23 0];
@btime gauss_jordan(A)

A = Float64[1 4 8  1
            4 5 6  8
            1 3 2  4
            4 5 67 23];
gauss_jordan(A)

# when dimensions doesn't match
A = [1 4 8 1  4;
     4 5 6 8 11;
     1 3 2 4  8];
gauss_jordan(A)


text:
function default_builtins()
    a = ""
    b = false
    c = 0
    @assert(a == "")
    @assert(b == false)
    @assert(c == 0)
end

function main()
    a = max(1, 2)
    println(join([a], " "))
    b = min(1, 2)
    println(join([b], " "))
end

main()


text:
##### version 2 -> including the leg

#=
renormalize!
=#

function renormalize!{T}(lattice::Classical2dFractalLattice{T}, dimM::Int)
	tensorT, tensorTtilde = getTensorT(lattice)

	tensorS1, tensorS1tilde, tensorSc1 =
	constructFirstHalf(tensorT, tensorTtilde, lattice, dimM, 2)

	constructSecondHalf(tensorS1, tensorS1tilde, lattice, dimM, 1, tensorSc1)
end

function calculateNextCoreTensor{T}(tensorT::Array{T,4}, tensorTtilde::Array{T,4}, dimM::Int; mirror = false)
	tenmatM = getTensorM(tensorT, mirror = mirror)
	matrixMMd = tenmatM.matrix * transpose(tenmatM.matrix)
	Ul, lambdaVector, Uld = svd(matrixMMd)
	trunUl = truncMatrixU(Ul,dimM)
	tensorUl = matU2tenU(trunUl)
	newTensorT = getNewTensorT(tensorT,tensorUl)
	newTensorTtilde = getNewTensorT(tensorT,tensorUl,tensorTtilde)
	return newTensorT, newTensorTtilde, tensorUl
end

function calculateScTensor{T}(tensorT::Array{T,4}, tensorCorner::Array{T,3}, tensorU::Array{T,3})
	@tensor begin
		dummy1[x1,x2,x1p,x2p,y] := tensorCorner[x1,x1p,s] * tensorT[x2p,x2,y,s]
		#####? maybe tensorU[x2,x1,i]
		dummy2[x,x1p,x2p,y] := tensorU[x1,x2,x] * dummy1[x1,x2,x1p,x2p,y]
		tensorSc[x,xp,y] := dummy2[x,x1p,x2p,y] * tensorU[x1p,x2p,xp]
	end
	return tensorSc
end

function calculateLegCaretTensor{T}(lattice::Classical2dFractalLattice{T}, XorY::Int, tensorU::Array{T,3})
	tensorP = getTensorP(lattice, XorY)
	legTensor = getNewLegTensor(tensorP, tensorU)
	caretTensor = getNewCaretTensor(getTensorQ(lattice), tensorU, XorY)

	legextension = getLegextension(lattice)
	if legextension == 0
		legTensor = legTensor^getLegextension(lattice)
		caretTensor = ones(size(caretTensor))
	elseif legextension > 0
		extraleg = legTensor^(getLegextension(lattice) - 1)
		legTensor = legTensor * extraleg
		caretTensor = extraleg * caretTensor
	else
		error("legextension should be not negative integer")
	end
	return legTensor, caretTensor
end

function constructFirstHalf{T}(tensorT::Array{T,4}, tensorTtilde::Array{T,4}, lattice::Classical2dFractalLattice{T}, dimM::Int, XorY::Int)
	newTensorT, newTensorTtilde, tensorU = calculateNextCoreTensor(tensorT, tensorTtilde, dimM)
	legTensor, caretTensor = calculateLegCaretTensor(lattice, XorY, tensorU)

	@tensor begin # put legs to one side
		tensorS[x,xp,y,yp] := # legTensor * newTensorT
		  legTensor[x,i] * newTensorT[i,xp,y,yp]
		tensorStilde[x,xp,y,yp] := # legTensor * newTensorTtilde
		  legTensor[x,i] * newTensorTtilde[i,xp,y,yp]
	  	tensorSc[y,yp,xp] := # caretTensor * newTensorT
			caretTensor[j] * newTensorT[j,xp,y,yp]
	end

	newRotTensorS = rotateTensorT90(tensorS)
	newRotTensorStilde = rotateTensorT90(tensorStilde)

	return newRotTensorS, newRotTensorStilde, tensorSc
end

function constructSecondHalf{T}(tensorT::Array{T,4}, tensorTtilde::Array{T,4}, lattice::Classical2dFractalLattice{T}, dimM::Int, XorY::Int, tensorCorner::Array{T,3})
	newTensorS, newTensorStilde, tensorU = calculateNextCoreTensor(tensorT, tensorTtilde, dimM, mirror = true)
	legTensor, caretTensor = calculateLegCaretTensor(lattice, XorY, tensorU)
	tensorSc = calculateScTensor(tensorT, tensorCorner, tensorU) # x, xp, y

	@tensor begin
		dummy1[x,xp,y,yp] := legTensor[x,i] * newTensorS[i,xp,y,yp] # legTensor * newTensorS
		newTensorT[x,xp,y,yp] := dummy1[x,j,y,yp] * legTensor[j,xp] # legTensor * newTensorS * legTensor = dummy1 * legTensor
		newTensorPy[y,yp,x] := dummy1[x,j,y,yp] * caretTensor[j] # legTensor * newTensorS * caretTensor = dummy1 * caretTensor
		dummy2[x,xp,yp] := legTensor[x,i] * tensorSc[i,xp,yp] # legTensor * ...
		newTensorPx[y,yp,s] := dummy2[y,j,s] * legTensor[j,yp]# legTensor * ... * leg
		newTensorQ[x,y] := dummy2[y,j,x] * caretTensor[j] # leg * ... * caret
		# tensorTtilde
		dummy1tilde[x,xp,y,yp] := legTensor[x,i] * newTensorStilde[i,xp,y,yp]
		newTensorTtilde[x,xp,y,yp] := dummy1tilde[x,j,y,yp] * legTensor[j,xp] # legTensor * newTensorS * legTensor = dummy1 * legTensor
	end

	newRotTensorT = rotateTensorT90(newTensorT)
	newRotTensorTtilde = rotateTensorT90(newTensorTtilde)

	setTensorT!(lattice, newRotTensorT)
	setTensorT!(lattice, newRotTensorTtilde; tilde = true)
	setTensorP!(lattice, 1, newTensorPx)
	setTensorP!(lattice, 2, newTensorPy)
	setTensorQ!(lattice, newTensorQ)

	return newRotTensorT, newRotTensorTtilde
end

function getNewLegTensor{T}(tensorP::Array{T,3}, tensorU::Array{T,3})
	@tensor begin
		#####? maybe tensorU[x2,x1,i]
		dummy1[i,x2,s,x1p] := tensorU[x1,x2,i] * tensorP[x1,x1p,s]
		dummy2[i,x1p,x2p] := dummy1[i,x2,s,x1p] * tensorP[x2p,x2,s]
		legTensor[i,j] :=  dummy2[i,x1p,x2p] * tensorU[x1p,x2p,j]
	end
	return legTensor
end

function getNewCaretTensor{T}(tensorQ::Array{T,2}, tensorU::Array{T,3}, XorY::Int)
	if XorY == 1 # C[X]
		@tensor begin
			dummy[y1,y2] := tensorQ[s,y1] * tensorQ[s,y2]
			caretTensor[j] := dummy[y1,y2] * tensorU[y1,y2,j]
		end
	else # C[Y]
		@tensor begin
			dummy[x1,x2] := tensorQ[x1,s] * tensorQ[x2,s]
			caretTensor[j] := dummy[x1,x2] * tensorU[x1,x2,j]
		end
	end
	return caretTensor
end


text:
#-----------------------------------------------------------------------# MLData
# something or nothing
const _N{T} = Union{Nothing, T} where T

struct MLData{X < _N{AbstractArray}, Y < _N{AbstractArray}, W <: _N{AbstractVector}}
    x::X
    y::Y
    w::W
    function MLData(x::X, y::Y, w::W) where {X,Y,W}
        new{X,Y,W}(x, y, w)
    end
end
# TODO: check that data sizes are correct in constructor
mldata(;x=nothing, y=nothing, w=nothing) = MLData(x, y, w)

function Base.show(io::IO, o::MLData)
    print(io, typeof(o))
    o.x != nothing && print(io, "\n    - x: ", summary(o.x))
    o.y != nothing && print(io, "\n    - y: ", summary(o.y))
    o.w != nothing && print(io, "\n    - w: ", summary(o.w))
end

xy(o::MLData) = o.x, o.y

nobs(o::MLData{<:AbstractArray}) = size(o.x, 1)
nobs(o::MLData{Nothing, <:AbstractArray}) = size(o,y, 1)

npredictors(o::MLData{<:AbstractArray}) = size(o.x, 2)
npredictors(o::MLData{Nothing}) = 0


text:
using Events
using Test

@testset "simple" begin
    x = 42
    ee = Emitter()
    off1! = on!(ee, :inc) do v
        x += v
    end
    off2! = on!(ee, :inc) do v
        x += v
    end
    @test x == 42
    emit(ee, :inc, 1)
    @test x == 42 + 2
    off2!()
    emit(ee, :inc, 1)
    @test x == 42 + 3
end

@testset "mixed types" begin
    x = 0
    ee = Emitter()
    f(v::Int) = x += v
    f(v::AbstractString) = x += parse(Int, v)
    on!(f, ee, :inc)
    emit(ee, :inc, 1)
    emit(ee, :inc, "1")
    @test x == 2
end


text:
int_rules_1_1_3_5 = @theory begin

    #= ::Subsection::Closed:: =#
    #= 1.1.3.5*(a+b*x^n)^p*(c+d*x^n)^q*(e+f*x^n)^r =#
    @apply_utils Antiderivative((~a + ~(b') * (~x) ^ ~n) ^ ~(p') * (~c + ~(d') * (~x) ^ ~n) ^ ~(q') * (~e + ~(f') * (~x) ^ ~n) ^ ~(r'), ~x) => Antiderivative(ExpandIntegrand((~a + ~b * (~x) ^ ~n) ^ ~p * (~c + ~d * (~x) ^ ~n) ^ ~q * (~e + ~f * (~x) ^ ~n) ^ ~r, ~x), ~x) <-- FreeQ([~a, ~b, ~c, ~d, ~e, ~f, ~n], ~x) && (IGtQ(~p, 0) && (IGtQ(~q, 0) && IGtQ(~r, 0)))
    @apply_utils Antiderivative((~e + ~(f') * (~x) ^ ~n) / ((~a + ~(b') * (~x) ^ ~n) * (~c + ~(d') * (~x) ^ ~n)), ~x) => ((~b * ~e - ~a * ~f) / (~b * ~c - ~a * ~d)) * Antiderivative(1 / (~a + ~b * (~x) ^ ~n), ~x) - ((~d * ~e - ~c * ~f) / (~b * ~c - ~a * ~d)) * Antiderivative(1 / (~c + ~d * (~x) ^ ~n), ~x) <-- FreeQ([~a, ~b, ~c, ~d, ~e, ~f, ~n], ~x)
    @apply_utils Antiderivative((~e + ~(f') * (~x) ^ ~n) / ((~a + ~(b') * (~x) ^ ~n) * sqrt(~c + ~(d') * (~x) ^ ~n)), ~x) => (~f / ~b) * Antiderivative(1 / sqrt(~c + ~d * (~x) ^ ~n), ~x) + ((~b * ~e - ~a * ~f) / ~b) * Antiderivative(1 / ((~a + ~b * (~x) ^ ~n) * sqrt(~c + ~d * (~x) ^ ~n)), ~x) <-- FreeQ([~a, ~b, ~c, ~d, ~e, ~f, ~n], ~x)
    @apply_utils Antiderivative((~e + ~(f') * (~x) ^ ~n) / (sqrt(~a + ~(b') * (~x) ^ ~n) * sqrt(~c + ~(d') * (~x) ^ ~n)), ~x) => (~f / ~b) * Antiderivative(sqrt(~a + ~b * (~x) ^ ~n) / sqrt(~c + ~d * (~x) ^ ~n), ~x) + ((~b * ~e - ~a * ~f) / ~b) * Antiderivative(1 / (sqrt(~a + ~b * (~x) ^ ~n) * sqrt(~c + ~d * (~x) ^ ~n)), ~x) <-- FreeQ([~a, ~b, ~c, ~d, ~e, ~f, ~n], ~x) && Not(EqQ(~n, 2) && (PosQ(~b / ~a) && PosQ(~d / ~c) || NegQ(~b / ~a) && (PosQ(~d / ~c) || GtQ(~a, 0) && (Not(GtQ(~c, 0)) || SimplerSqrtQ(-(~b) / ~a, -(~d) / ~c)))))
    @apply_utils Antiderivative((~e + ~(f') * (~x) ^ 2) / (sqrt(~a + ~(b') * (~x) ^ 2) * (~c + ~(d') * (~x) ^ 2) ^ (3 / 2)), ~x) => ((~b * ~e - ~a * ~f) / (~b * ~c - ~a * ~d)) * Antiderivative(1 / (sqrt(~a + ~b * (~x) ^ 2) * sqrt(~c + ~d * (~x) ^ 2)), ~x) - ((~d * ~e - ~c * ~f) / (~b * ~c - ~a * ~d)) * Antiderivative(sqrt(~a + ~b * (~x) ^ 2) / (~c + ~d * (~x) ^ 2) ^ (3 / 2), ~x) <-- FreeQ([~a, ~b, ~c, ~d, ~e, ~f], ~x) && (PosQ(~b / ~a) && PosQ(~d / ~c))
    @apply_utils Antiderivative((~a + ~(b') * (~x) ^ ~n) ^ ~p * (~c + ~(d') * (~x) ^ ~n) ^ ~(q') * (~e + ~(f') * (~x) ^ ~n), ~x) => (-((~b * ~e - ~a * ~f)) * ~x * (~a + ~b * (~x) ^ ~n) ^ (~p + 1) * (~c + ~d * (~x) ^ ~n) ^ ~q) / (~a * ~b * ~n * (~p + 1)) + (1 / (~a * ~b * ~n * (~p + 1))) * Antiderivative((~a + ~b * (~x) ^ ~n) ^ (~p + 1) * (~c + ~d * (~x) ^ ~n) ^ (~q - 1) * Simp(~c * ((~b * ~e * ~n * (~p + 1) + ~b * ~e) - ~a * ~f) + ~d * (~b * ~e * ~n * (~p + 1) + (~b * ~e - ~a * ~f) * (~n * ~q + 1)) * (~x) ^ ~n, ~x), ~x) <-- FreeQ([~a, ~b, ~c, ~d, ~e, ~f, ~n], ~x) && (LtQ(~p, -1) && GtQ(~q, 0))
    @apply_utils Antiderivative((~a + ~(b') * (~x) ^ ~n) ^ ~p * (~c + ~(d') * (~x) ^ ~n) ^ ~(q') * (~e + ~(f') * (~x) ^ ~n), ~x) => (-((~b * ~e - ~a * ~f)) * ~x * (~a + ~b * (~x) ^ ~n) ^ (~p + 1) * (~c + ~d * (~x) ^ ~n) ^ (~q + 1)) / (~a * ~n * (~b * ~c - ~a * ~d) * (~p + 1)) + (1 / (~a * ~n * (~b * ~c - ~a * ~d) * (~p + 1))) * Antiderivative((~a + ~b * (~x) ^ ~n) ^ (~p + 1) * (~c + ~d * (~x) ^ ~n) ^ ~q * Simp(~c * (~b * ~e - ~a * ~f) + ~e * ~n * (~b * ~c - ~a * ~d) * (~p + 1) + ~d * (~b * ~e - ~a * ~f) * (~n * (~p + ~q + 2) + 1) * (~x) ^ ~n, ~x), ~x) <-- FreeQ([~a, ~b, ~c, ~d, ~e, ~f, ~n, ~q], ~x) && LtQ(~p, -1)
    @apply_utils Antiderivative((~a + ~(b') * (~x) ^ ~n) ^ ~(p') * (~c + ~(d') * (~x) ^ ~n) ^ ~(q') * (~e + ~(f') * (~x) ^ ~n), ~x) => (~f * ~x * (~a + ~b * (~x) ^ ~n) ^ (~p + 1) * (~c + ~d * (~x) ^ ~n) ^ ~q) / (~b * (~n * (~p + ~q + 1) + 1)) + (1 / (~b * (~n * (~p + ~q + 1) + 1))) * Antiderivative((~a + ~b * (~x) ^ ~n) ^ ~p * (~c + ~d * (~x) ^ ~n) ^ (~q - 1) * Simp(~c * ((~b * ~e - ~a * ~f) + ~b * ~e * ~n * (~p + ~q + 1)) + (~d * (~b * ~e - ~a * ~f) + ~f * ~n * ~q * (~b * ~c - ~a * ~d) + ~b * ~d * ~e * ~n * (~p + ~q + 1)) * (~x) ^ ~n, ~x), ~x) <-- FreeQ([~a, ~b, ~c, ~d, ~e, ~f, ~n, ~p], ~x) && (GtQ(~q, 0) && NeQ(~n * (~p + ~q + 1) + 1, 0))
    @apply_utils Antiderivative((~e + ~(f') * (~x) ^ 4) / ((~a + ~(b') * (~x) ^ 4) ^ (3 / 4) * (~c + ~(d') * (~x) ^ 4)), ~x) => ((~b * ~e - ~a * ~f) / (~b * ~c - ~a * ~d)) * Antiderivative(1 / (~a + ~b * (~x) ^ 4) ^ (3 / 4), ~x) - ((~d * ~e - ~c * ~f) / (~b * ~c - ~a * ~d)) * Antiderivative((~a + ~b * (~x) ^ 4) ^ (1 / 4) / (~c + ~d * (~x) ^ 4), ~x) <-- FreeQ([~a, ~b, ~c, ~d, ~e, ~f], ~x)
    @apply_utils Antiderivative(((~a + ~(b') * (~x) ^ ~n) ^ ~p * (~e + ~(f') * (~x) ^ ~n)) / (~c + ~(d') * (~x) ^ ~n), ~x) => (~f / ~d) * Antiderivative((~a + ~b * (~x) ^ ~n) ^ ~p, ~x) + ((~d * ~e - ~c * ~f) / ~d) * Antiderivative((~a + ~b * (~x) ^ ~n) ^ ~p / (~c + ~d * (~x) ^ ~n), ~x) <-- FreeQ([~a, ~b, ~c, ~d, ~e, ~f, ~p, ~n], ~x)
    @apply_utils Antiderivative((~a + ~(b') * (~x) ^ ~n) ^ ~(p') * (~c + ~(d') * (~x) ^ ~n) ^ ~(q') * (~e + ~(f') * (~x) ^ ~n), ~x) => ~e * Antiderivative((~a + ~b * (~x) ^ ~n) ^ ~p * (~c + ~d * (~x) ^ ~n) ^ ~q, ~x) + ~f * Antiderivative((~x) ^ ~n * (~a + ~b * (~x) ^ ~n) ^ ~p * (~c + ~d * (~x) ^ ~n) ^ ~q, ~x) <-- FreeQ([~a, ~b, ~c, ~d, ~e, ~f, ~n, ~p, ~q], ~x)
    @apply_utils Antiderivative(1 / ((~a + ~(b') * (~x) ^ 2) * (~c + ~(d') * (~x) ^ 2) * sqrt(~e + ~(f') * (~x) ^ 2)), ~x) => (~b / (~b * ~c - ~a * ~d)) * Antiderivative(1 / ((~a + ~b * (~x) ^ 2) * sqrt(~e + ~f * (~x) ^ 2)), ~x) - (~d / (~b * ~c - ~a * ~d)) * Antiderivative(1 / ((~c + ~d * (~x) ^ 2) * sqrt(~e + ~f * (~x) ^ 2)), ~x) <-- FreeQ([~a, ~b, ~c, ~d, ~e, ~f], ~x)
    @apply_utils Antiderivative(1 / ((~x) ^ 2 * (~c + ~(d') * (~x) ^ 2) * sqrt(~e + ~(f') * (~x) ^ 2)), ~x) => (1 / ~c) * Antiderivative(1 / ((~x) ^ 2 * sqrt(~e + ~f * (~x) ^ 2)), ~x) - (~d / ~c) * Antiderivative(1 / ((~c + ~d * (~x) ^ 2) * sqrt(~e + ~f * (~x) ^ 2)), ~x) <-- FreeQ([~c, ~d, ~e, ~f], ~x) && NeQ(~d * ~e - ~c * ~f, 0)
    @apply_utils Antiderivative((sqrt(~c + ~(d') * (~x) ^ 2) * sqrt(~e + ~(f') * (~x) ^ 2)) / (~a + ~(b') * (~x) ^ 2), ~x) => (~d / ~b) * Antiderivative(sqrt(~e + ~f * (~x) ^ 2) / sqrt(~c + ~d * (~x) ^ 2), ~x) + ((~b * ~c - ~a * ~d) / ~b) * Antiderivative(sqrt(~e + ~f * (~x) ^ 2) / ((~a + ~b * (~x) ^ 2) * sqrt(~c + ~d * (~x) ^ 2)), ~x) <-- FreeQ([~a, ~b, ~c, ~d, ~e, ~f], ~x) && (GtQ(~d / ~c, 0) && (GtQ(~f / ~e, 0) && Not(SimplerSqrtQ(~d / ~c, ~f / ~e))))
    @apply_utils Antiderivative((sqrt(~c + ~(d') * (~x) ^ 2) * sqrt(~e + ~(f') * (~x) ^ 2)) / (~a + ~(b') * (~x) ^ 2), ~x) => (~d / ~b) * Antiderivative(sqrt(~e + ~f * (~x) ^ 2) / sqrt(~c + ~d * (~x) ^ 2), ~x) + ((~b * ~c - ~a * ~d) / ~b) * Antiderivative(sqrt(~e + ~f * (~x) ^ 2) / ((~a + ~b * (~x) ^ 2) * sqrt(~c + ~d * (~x) ^ 2)), ~x) <-- FreeQ([~a, ~b, ~c, ~d, ~e, ~f], ~x) && Not(SimplerSqrtQ(-(~f) / ~e, -(~d) / ~c))
    @apply_utils Antiderivative(1 / ((~a + ~(b') * (~x) ^ 2) * sqrt(~c + ~(d') * (~x) ^ 2) * sqrt(~e + ~(f') * (~x) ^ 2)), ~x) => (-(~f) / (~b * ~e - ~a * ~f)) * Antiderivative(1 / (sqrt(~c + ~d * (~x) ^ 2) * sqrt(~e + ~f * (~x) ^ 2)), ~x) + (~b / (~b * ~e - ~a * ~f)) * Antiderivative(sqrt(~e + ~f * (~x) ^ 2) / ((~a + ~b * (~x) ^ 2) * sqrt(~c + ~d * (~x) ^ 2)), ~x) <-- FreeQ([~a, ~b, ~c, ~d, ~e, ~f], ~x) && (GtQ(~d / ~c, 0) && (GtQ(~f / ~e, 0) && Not(SimplerSqrtQ(~d / ~c, ~f / ~e))))
    @apply_utils Antiderivative(1 / ((~a + ~(b') * (~x) ^ 2) * sqrt(~c + ~(d') * (~x) ^ 2) * sqrt(~e + ~(f') * (~x) ^ 2)), ~x) => (1 / (~a * sqrt(~c) * sqrt(~e) * Rt(-(~d) / ~c, 2))) * Elliptic.Pi((~b * ~c) / (~a * ~d), asin(Rt(-(~d) / ~c, 2) * ~x), (~c * ~f) / (~d * ~e)) <-- FreeQ([~a, ~b, ~c, ~d, ~e, ~f], ~x) && (Not(GtQ(~d / ~c, 0)) && (GtQ(~c, 0) && (GtQ(~e, 0) && Not(Not(GtQ(~f / ~e, 0)) && SimplerSqrtQ(-(~f) / ~e, -(~d) / ~c)))))
    @apply_utils Antiderivative(1 / ((~a + ~(b') * (~x) ^ 2) * sqrt(~c + ~(d') * (~x) ^ 2) * sqrt(~e + ~(f') * (~x) ^ 2)), ~x) => (sqrt(1 + (~d / ~c) * (~x) ^ 2) / sqrt(~c + ~d * (~x) ^ 2)) * Antiderivative(1 / ((~a + ~b * (~x) ^ 2) * sqrt(1 + (~d / ~c) * (~x) ^ 2) * sqrt(~e + ~f * (~x) ^ 2)), ~x) <-- FreeQ([~a, ~b, ~c, ~d, ~e, ~f], ~x) && Not(GtQ(~c, 0))
    @apply_utils Antiderivative(sqrt(~c + ~(d') * (~x) ^ 2) / ((~a + ~(b') * (~x) ^ 2) * sqrt(~e + ~(f') * (~x) ^ 2)), ~x) => ((~c * sqrt(~e + ~f * (~x) ^ 2)) / (~a * ~e * Rt(~d / ~c, 2) * sqrt(~c + ~d * (~x) ^ 2) * sqrt((~c * (~e + ~f * (~x) ^ 2)) / (~e * (~c + ~d * (~x) ^ 2))))) * Elliptic.Pi(1 - (~b * ~c) / (~a * ~d), atan(Rt(~d / ~c, 2) * ~x), 1 - (~c * ~f) / (~d * ~e)) <-- FreeQ([~a, ~b, ~c, ~d, ~e, ~f], ~x) && PosQ(~d / ~c)
    #= Antiderivative(sqrt((~c)+(~d')*(~x)^2)/(((~a)+(~b')*(~x)^2)*sqrt((~e)+(~f')*(~x)^2)),~x) :=   sqrt(c+d*x^2)*sqrt(c*(e+f*x^2)/(e*(c+d*x^2)))/(a*Rt(d/c,2)*sqrt(e+f* x^2))* Elliptic.Pi(1-b*c/(a*d),atan(Rt(d/c,2)*x),1-c*f/(d*e)) <--  FreeQ([a,b,c,d,e,f],x) && PosQ(d/c) =#
    @apply_utils Antiderivative(sqrt(~c + ~(d') * (~x) ^ 2) / ((~a + ~(b') * (~x) ^ 2) * sqrt(~e + ~(f') * (~x) ^ 2)), ~x) => (~d / ~b) * Antiderivative(1 / (sqrt(~c + ~d * (~x) ^ 2) * sqrt(~e + ~f * (~x) ^ 2)), ~x) + ((~b * ~c - ~a * ~d) / ~b) * Antiderivative(1 / ((~a + ~b * (~x) ^ 2) * sqrt(~c + ~d * (~x) ^ 2) * sqrt(~e + ~f * (~x) ^ 2)), ~x) <-- FreeQ([~a, ~b, ~c, ~d, ~e, ~f], ~x) && NegQ(~d / ~c)
    @apply_utils Antiderivative(sqrt(~e + ~(f') * (~x) ^ 2) / ((~a + ~(b') * (~x) ^ 2) * (~c + ~(d') * (~x) ^ 2) ^ (3 / 2)), ~x) => (~b / (~b * ~c - ~a * ~d)) * Antiderivative(sqrt(~e + ~f * (~x) ^ 2) / ((~a + ~b * (~x) ^ 2) * sqrt(~c + ~d * (~x) ^ 2)), ~x) - (~d / (~b * ~c - ~a * ~d)) * Antiderivative(sqrt(~e + ~f * (~x) ^ 2) / (~c + ~d * (~x) ^ 2) ^ (3 / 2), ~x) <-- FreeQ([~a, ~b, ~c, ~d, ~e, ~f], ~x) && (PosQ(~d / ~c) && PosQ(~f / ~e))
    @apply_utils Antiderivative((~e + ~(f') * (~x) ^ 2) ^ (3 / 2) / ((~a + ~(b') * (~x) ^ 2) * (~c + ~(d') * (~x) ^ 2) ^ (3 / 2)), ~x) => ((~b * ~e - ~a * ~f) / (~b * ~c - ~a * ~d)) * Antiderivative(sqrt(~e + ~f * (~x) ^ 2) / ((~a + ~b * (~x) ^ 2) * sqrt(~c + ~d * (~x) ^ 2)), ~x) - ((~d * ~e - ~c * ~f) / (~b * ~c - ~a * ~d)) * Antiderivative(sqrt(~e + ~f * (~x) ^ 2) / (~c + ~d * (~x) ^ 2) ^ (3 / 2), ~x) <-- FreeQ([~a, ~b, ~c, ~d, ~e, ~f], ~x) && (PosQ(~d / ~c) && PosQ(~f / ~e))
    @apply_utils Antiderivative(((~c + ~(d') * (~x) ^ 2) ^ (3 / 2) * sqrt(~e + ~(f') * (~x) ^ 2)) / (~a + ~(b') * (~x) ^ 2), ~x) => ((~b * ~c - ~a * ~d) ^ 2 / (~b) ^ 2) * Antiderivative(sqrt(~e + ~f * (~x) ^ 2) / ((~a + ~b * (~x) ^ 2) * sqrt(~c + ~d * (~x) ^ 2)), ~x) + (~d / (~b) ^ 2) * Antiderivative((((2 * ~b * ~c - ~a * ~d) + ~b * ~d * (~x) ^ 2) * sqrt(~e + ~f * (~x) ^ 2)) / sqrt(~c + ~d * (~x) ^ 2), ~x) <-- FreeQ([~a, ~b, ~c, ~d, ~e, ~f], ~x) && (PosQ(~d / ~c) && PosQ(~f / ~e))
    @apply_utils Antiderivative(((~c + ~(d') * (~x) ^ 2) ^ ~q * (~e + ~(f') * (~x) ^ 2) ^ ~r) / (~a + ~(b') * (~x) ^ 2), ~x) => ((~b * (~b * ~e - ~a * ~f)) / (~b * ~c - ~a * ~d) ^ 2) * Antiderivative(((~c + ~d * (~x) ^ 2) ^ (~q + 2) * (~e + ~f * (~x) ^ 2) ^ (~r - 1)) / (~a + ~b * (~x) ^ 2), ~x) - (1 / (~b * ~c - ~a * ~d) ^ 2) * Antiderivative((~c + ~d * (~x) ^ 2) ^ ~q * (~e + ~f * (~x) ^ 2) ^ (~r - 1) * (((2 * ~b * ~c * ~d * ~e - ~a * (~d) ^ 2 * ~e) - ~b * (~c) ^ 2 * ~f) + (~d) ^ 2 * (~b * ~e - ~a * ~f) * (~x) ^ 2), ~x) <-- FreeQ([~a, ~b, ~c, ~d, ~e, ~f], ~x) && (LtQ(~q, -1) && GtQ(~r, 1))
    @apply_utils Antiderivative(((~c + ~(d') * (~x) ^ 2) ^ ~q * (~e + ~(f') * (~x) ^ 2) ^ ~r) / (~a + ~(b') * (~x) ^ 2), ~x) => (~d / ~b) * Antiderivative((~c + ~d * (~x) ^ 2) ^ (~q - 1) * (~e + ~f * (~x) ^ 2) ^ ~r, ~x) + ((~b * ~c - ~a * ~d) / ~b) * Antiderivative(((~c + ~d * (~x) ^ 2) ^ (~q - 1) * (~e + ~f * (~x) ^ 2) ^ ~r) / (~a + ~b * (~x) ^ 2), ~x) <-- FreeQ([~a, ~b, ~c, ~d, ~e, ~f, ~r], ~x) && GtQ(~q, 1)
    @apply_utils Antiderivative(((~c + ~(d') * (~x) ^ 2) ^ ~q * (~e + ~(f') * (~x) ^ 2) ^ ~r) / (~a + ~(b') * (~x) ^ 2), ~x) => ((~b) ^ 2 / (~b * ~c - ~a * ~d) ^ 2) * Antiderivative(((~c + ~d * (~x) ^ 2) ^ (~q + 2) * (~e + ~f * (~x) ^ 2) ^ ~r) / (~a + ~b * (~x) ^ 2), ~x) - (~d / (~b * ~c - ~a * ~d) ^ 2) * Antiderivative((~c + ~d * (~x) ^ 2) ^ ~q * (~e + ~f * (~x) ^ 2) ^ ~r * ((2 * ~b * ~c - ~a * ~d) + ~b * ~d * (~x) ^ 2), ~x) <-- FreeQ([~a, ~b, ~c, ~d, ~e, ~f, ~r], ~x) && LtQ(~q, -1)
    @apply_utils Antiderivative(((~c + ~(d') * (~x) ^ 2) ^ ~q * (~e + ~(f') * (~x) ^ 2) ^ ~r) / (~a + ~(b') * (~x) ^ 2), ~x) => (-(~d) / (~b * ~c - ~a * ~d)) * Antiderivative((~c + ~d * (~x) ^ 2) ^ ~q * (~e + ~f * (~x) ^ 2) ^ ~r, ~x) + (~b / (~b * ~c - ~a * ~d)) * Antiderivative(((~c + ~d * (~x) ^ 2) ^ (~q + 1) * (~e + ~f * (~x) ^ 2) ^ ~r) / (~a + ~b * (~x) ^ 2), ~x) <-- FreeQ([~a, ~b, ~c, ~d, ~e, ~f, ~r], ~x) && LeQ(~q, -1)
    @apply_utils Antiderivative((sqrt(~c + ~(d') * (~x) ^ 2) * sqrt(~e + ~(f') * (~x) ^ 2)) / (~a + ~(b') * (~x) ^ 2) ^ 2, ~x) => (~x * sqrt(~c + ~d * (~x) ^ 2) * sqrt(~e + ~f * (~x) ^ 2)) / (2 * ~a * (~a + ~b * (~x) ^ 2)) + ((~d * ~f) / (2 * ~a * (~b) ^ 2)) * Antiderivative((~a - ~b * (~x) ^ 2) / (sqrt(~c + ~d * (~x) ^ 2) * sqrt(~e + ~f * (~x) ^ 2)), ~x) + (((~b) ^ 2 * ~c * ~e - (~a) ^ 2 * ~d * ~f) / (2 * ~a * (~b) ^ 2)) * Antiderivative(1 / ((~a + ~b * (~x) ^ 2) * sqrt(~c + ~d * (~x) ^ 2) * sqrt(~e + ~f * (~x) ^ 2)), ~x) <-- FreeQ([~a, ~b, ~c, ~d, ~e, ~f], ~x)
    @apply_utils Antiderivative(1 / ((~a + ~(b') * (~x) ^ 2) ^ 2 * sqrt(~c + ~(d') * (~x) ^ 2) * sqrt(~e + ~(f') * (~x) ^ 2)), ~x) => (((~b) ^ 2 * ~x * sqrt(~c + ~d * (~x) ^ 2) * sqrt(~e + ~f * (~x) ^ 2)) / (2 * ~a * (~b * ~c - ~a * ~d) * (~b * ~e - ~a * ~f) * (~a + ~b * (~x) ^ 2)) - ((~d * ~f) / (2 * ~a * (~b * ~c - ~a * ~d) * (~b * ~e - ~a * ~f))) * Antiderivative((~a + ~b * (~x) ^ 2) / (sqrt(~c + ~d * (~x) ^ 2) * sqrt(~e + ~f * (~x) ^ 2)), ~x)) + ((((~b) ^ 2 * ~c * ~e + 3 * (~a) ^ 2 * ~d * ~f) - 2 * ~a * ~b * (~d * ~e + ~c * ~f)) / (2 * ~a * (~b * ~c - ~a * ~d) * (~b * ~e - ~a * ~f))) * Antiderivative(1 / ((~a + ~b * (~x) ^ 2) * sqrt(~c + ~d * (~x) ^ 2) * sqrt(~e + ~f * (~x) ^ 2)), ~x) <-- FreeQ([~a, ~b, ~c, ~d, ~e, ~f], ~x)
    @apply_utils Antiderivative((~a + ~(b') * (~x) ^ ~n) ^ ~p * (~c + ~(d') * (~x) ^ ~n) ^ ~q * (~e + ~(f') * (~x) ^ ~n) ^ ~r, ~x) => (~d / ~b) * Antiderivative((~a + ~b * (~x) ^ ~n) ^ (~p + 1) * (~c + ~d * (~x) ^ ~n) ^ (~q - 1) * (~e + ~f * (~x) ^ ~n) ^ ~r, ~x) + ((~b * ~c - ~a * ~d) / ~b) * Antiderivative((~a + ~b * (~x) ^ ~n) ^ ~p * (~c + ~d * (~x) ^ ~n) ^ (~q - 1) * (~e + ~f * (~x) ^ ~n) ^ ~r, ~x) <-- FreeQ([~a, ~b, ~c, ~d, ~e, ~f, ~n, ~r], ~x) && (ILtQ(~p, 0) && GtQ(~q, 0))
    @apply_utils Antiderivative((~a + ~(b') * (~x) ^ ~n) ^ ~p * (~c + ~(d') * (~x) ^ ~n) ^ ~q * (~e + ~(f') * (~x) ^ ~n) ^ ~r, ~x) => (~b / (~b * ~c - ~a * ~d)) * Antiderivative((~a + ~b * (~x) ^ ~n) ^ ~p * (~c + ~d * (~x) ^ ~n) ^ (~q + 1) * (~e + ~f * (~x) ^ ~n) ^ ~r, ~x) - (~d / (~b * ~c - ~a * ~d)) * Antiderivative((~a + ~b * (~x) ^ ~n) ^ (~p + 1) * (~c + ~d * (~x) ^ ~n) ^ ~q * (~e + ~f * (~x) ^ ~n) ^ ~r, ~x) <-- FreeQ([~a, ~b, ~c, ~d, ~e, ~f, ~n, ~q], ~x) && (ILtQ(~p, 0) && LeQ(~q, -1))
    @apply_utils Antiderivative(1 / (sqrt(~a + ~(b') * (~x) ^ 2) * sqrt(~c + ~(d') * (~x) ^ 2) * sqrt(~e + ~(f') * (~x) ^ 2)), ~x) => ((sqrt(~c + ~d * (~x) ^ 2) * sqrt((~a * (~e + ~f * (~x) ^ 2)) / (~e * (~a + ~b * (~x) ^ 2)))) / (~c * sqrt(~e + ~f * (~x) ^ 2) * sqrt((~a * (~c + ~d * (~x) ^ 2)) / (~c * (~a + ~b * (~x) ^ 2))))) * Subst(Antiderivative(1 / (sqrt(1 - ((~b * ~c - ~a * ~d) * (~x) ^ 2) / ~c) * sqrt(1 - ((~b * ~e - ~a * ~f) * (~x) ^ 2) / ~e)), ~x), ~x, ~x / sqrt(~a + ~b * (~x) ^ 2)) <-- FreeQ([~a, ~b, ~c, ~d, ~e, ~f], ~x)
    @apply_utils Antiderivative(sqrt(~a + ~(b') * (~x) ^ 2) / (sqrt(~c + ~(d') * (~x) ^ 2) * sqrt(~e + ~(f') * (~x) ^ 2)), ~x) => ((~a * sqrt(~c + ~d * (~x) ^ 2) * sqrt((~a * (~e + ~f * (~x) ^ 2)) / (~e * (~a + ~b * (~x) ^ 2)))) / (~c * sqrt(~e + ~f * (~x) ^ 2) * sqrt((~a * (~c + ~d * (~x) ^ 2)) / (~c * (~a + ~b * (~x) ^ 2))))) * Subst(Antiderivative(1 / ((1 - ~b * (~x) ^ 2) * sqrt(1 - ((~b * ~c - ~a * ~d) * (~x) ^ 2) / ~c) * sqrt(1 - ((~b * ~e - ~a * ~f) * (~x) ^ 2) / ~e)), ~x), ~x, ~x / sqrt(~a + ~b * (~x) ^ 2)) <-- FreeQ([~a, ~b, ~c, ~d, ~e, ~f], ~x)
    @apply_utils Antiderivative(sqrt(~c + ~(d') * (~x) ^ 2) / ((~a + ~(b') * (~x) ^ 2) ^ (3 / 2) * sqrt(~e + ~(f') * (~x) ^ 2)), ~x) => ((sqrt(~c + ~d * (~x) ^ 2) * sqrt((~a * (~e + ~f * (~x) ^ 2)) / (~e * (~a + ~b * (~x) ^ 2)))) / (~a * sqrt(~e + ~f * (~x) ^ 2) * sqrt((~a * (~c + ~d * (~x) ^ 2)) / (~c * (~a + ~b * (~x) ^ 2))))) * Subst(Antiderivative(sqrt(1 - ((~b * ~c - ~a * ~d) * (~x) ^ 2) / ~c) / sqrt(1 - ((~b * ~e - ~a * ~f) * (~x) ^ 2) / ~e), ~x), ~x, ~x / sqrt(~a + ~b * (~x) ^ 2)) <-- FreeQ([~a, ~b, ~c, ~d, ~e, ~f], ~x)
    @apply_utils Antiderivative((sqrt(~a + ~(b') * (~x) ^ 2) * sqrt(~c + ~(d') * (~x) ^ 2)) / sqrt(~e + ~(f') * (~x) ^ 2), ~x) => (((~d * ~x * sqrt(~a + ~b * (~x) ^ 2) * sqrt(~e + ~f * (~x) ^ 2)) / (2 * ~f * sqrt(~c + ~d * (~x) ^ 2)) - ((~c * (~d * ~e - ~c * ~f)) / (2 * ~f)) * Antiderivative(sqrt(~a + ~b * (~x) ^ 2) / ((~c + ~d * (~x) ^ 2) ^ (3 / 2) * sqrt(~e + ~f * (~x) ^ 2)), ~x)) + ((~b * ~c * (~d * ~e - ~c * ~f)) / (2 * ~d * ~f)) * Antiderivative(1 / (sqrt(~a + ~b * (~x) ^ 2) * sqrt(~c + ~d * (~x) ^ 2) * sqrt(~e + ~f * (~x) ^ 2)), ~x)) - (((~b * ~d * ~e - ~b * ~c * ~f) - ~a * ~d * ~f) / (2 * ~d * ~f)) * Antiderivative(sqrt(~c + ~d * (~x) ^ 2) / (sqrt(~a + ~b * (~x) ^ 2) * sqrt(~e + ~f * (~x) ^ 2)), ~x) <-- FreeQ([~a, ~b, ~c, ~d, ~e, ~f], ~x) && PosQ((~d * ~e - ~c * ~f) / ~c)
    @apply_utils Antiderivative((sqrt(~a + ~(b') * (~x) ^ 2) * sqrt(~c + ~(d') * (~x) ^ 2)) / sqrt(~e + ~(f') * (~x) ^ 2), ~x) => ((~x * sqrt(~a + ~b * (~x) ^ 2) * sqrt(~c + ~d * (~x) ^ 2)) / (2 * sqrt(~e + ~f * (~x) ^ 2)) + ((~e * (~b * ~e - ~a * ~f)) / (2 * ~f)) * Antiderivative(sqrt(~c + ~d * (~x) ^ 2) / (sqrt(~a + ~b * (~x) ^ 2) * (~e + ~f * (~x) ^ 2) ^ (3 / 2)), ~x) + (((~b * ~e - ~a * ~f) * (~d * ~e - 2 * ~c * ~f)) / (2 * (~f) ^ 2)) * Antiderivative(1 / (sqrt(~a + ~b * (~x) ^ 2) * sqrt(~c + ~d * (~x) ^ 2) * sqrt(~e + ~f * (~x) ^ 2)), ~x)) - (((~b * ~d * ~e - ~b * ~c * ~f) - ~a * ~d * ~f) / (2 * (~f) ^ 2)) * Antiderivative(sqrt(~e + ~f * (~x) ^ 2) / (sqrt(~a + ~b * (~x) ^ 2) * sqrt(~c + ~d * (~x) ^ 2)), ~x) <-- FreeQ([~a, ~b, ~c, ~d, ~e, ~f], ~x) && NegQ((~d * ~e - ~c * ~f) / ~c)
    @apply_utils Antiderivative((sqrt(~a + ~(b') * (~x) ^ 2) * sqrt(~c + ~(d') * (~x) ^ 2)) / (~e + ~(f') * (~x) ^ 2) ^ (3 / 2), ~x) => (~b / ~f) * Antiderivative(sqrt(~c + ~d * (~x) ^ 2) / (sqrt(~a + ~b * (~x) ^ 2) * sqrt(~e + ~f * (~x) ^ 2)), ~x) - ((~b * ~e - ~a * ~f) / ~f) * Antiderivative(sqrt(~c + ~d * (~x) ^ 2) / (sqrt(~a + ~b * (~x) ^ 2) * (~e + ~f * (~x) ^ 2) ^ (3 / 2)), ~x) <-- FreeQ([~a, ~b, ~c, ~d, ~e, ~f], ~x)
    @apply_utils Antiderivative((~a + ~(b') * (~x) ^ ~n) ^ ~p * (~c + ~(d') * (~x) ^ ~n) ^ ~q * (~e + ~(f') * (~x) ^ ~n) ^ ~r, ~x) => With([u = ExpandIntegrand((~a + ~b * (~x) ^ ~n) ^ ~p * (~c + ~d * (~x) ^ ~n) ^ ~q * (~e + ~f * (~x) ^ ~n) ^ ~r, ~x)], Antiderivative(u, ~x) <-- SumQ(u)) <-- FreeQ([~a, ~b, ~c, ~d, ~e, ~f, ~p, ~q, ~r], ~x) && IGtQ(~n, 0)
    @apply_utils Antiderivative((~a + ~(b') * (~x) ^ ~n) ^ ~p * (~c + ~(d') * (~x) ^ ~n) ^ ~q * (~e + ~(f') * (~x) ^ ~n) ^ ~r, ~x) => -(Subst(Antiderivative(((~a + ~b * (~x) ^ -(~n)) ^ ~p * (~c + ~d * (~x) ^ -(~n)) ^ ~q * (~e + ~f * (~x) ^ -(~n)) ^ ~r) / (~x) ^ 2, ~x), ~x, 1 / ~x)) <-- FreeQ([~a, ~b, ~c, ~d, ~e, ~f, ~p, ~q, ~r], ~x) && ILtQ(~n, 0)
    @apply_utils Antiderivative((~a + ~(b') * (~x) ^ ~n) ^ ~(p') * (~c + ~(d') * (~x) ^ ~n) ^ ~(q') * (~e + ~(f') * (~x) ^ ~n) ^ ~(r'), ~x) => Unintegrable((~a + ~b * (~x) ^ ~n) ^ ~p * (~c + ~d * (~x) ^ ~n) ^ ~q * (~e + ~f * (~x) ^ ~n) ^ ~r, ~x) <-- FreeQ([~a, ~b, ~c, ~d, ~e, ~f, ~n, ~p, ~q, ~r], ~x)
    @apply_utils Antiderivative((~(a') + ~(b') * (~u) ^ ~n) ^ ~(p') * (~(c') + ~(d') * (~v) ^ ~n) ^ ~(q') * (~(e') + ~(f') * (~w) ^ ~n) ^ ~(r'), ~x) => (1 / Coefficient(~u, ~x, 1)) * Subst(Antiderivative((~a + ~b * (~x) ^ ~n) ^ ~p * (~c + ~d * (~x) ^ ~n) ^ ~q * (~e + ~f * (~x) ^ ~n) ^ ~r, ~x), ~x, ~u) <-- FreeQ([~a, ~b, ~c, ~d, ~e, ~f, ~p, ~n, ~q, ~r], ~x) && (EqQ(~u, ~v) && (EqQ(~u, ~w) && (LinearQ(~u, ~x) && NeQ(~u, ~x))))
    @apply_utils Antiderivative((~(a') + ~(b') * (~x) ^ ~(n')) ^ ~(p') * (~c + ~(d') * (~x) ^ ~(mn')) ^ ~(q') * (~e + ~(f') * (~x) ^ ~(n')) ^ ~(r'), ~x) => Antiderivative(((~a + ~b * (~x) ^ ~n) ^ ~p * (~d + ~c * (~x) ^ ~n) ^ ~q * (~e + ~f * (~x) ^ ~n) ^ ~r) / (~x) ^ (~n * ~q), ~x) <-- FreeQ([~a, ~b, ~c, ~d, ~e, ~f, ~n, ~p, ~r], ~x) && (EqQ(~mn, -(~n)) && IntegerQ(~q))
    @apply_utils Antiderivative((~(a') + ~(b') * (~x) ^ ~(n')) ^ ~(p') * (~c + ~(d') * (~x) ^ ~(mn')) ^ ~(q') * (~e + ~(f') * (~x) ^ ~(n')) ^ ~(r'), ~x) => Antiderivative((~x) ^ (~n * (~p + ~r)) * (~b + ~a * (~x) ^ -(~n)) ^ ~p * (~c + ~d * (~x) ^ -(~n)) ^ ~q * (~f + ~e * (~x) ^ -(~n)) ^ ~r, ~x) <-- FreeQ([~a, ~b, ~c, ~d, ~e, ~f, ~n, ~q], ~x) && (EqQ(~mn, -(~n)) && (IntegerQ(~p) && IntegerQ(~r)))
    @apply_utils Antiderivative((~(a') + ~(b') * (~x) ^ ~(n')) ^ ~(p') * (~c + ~(d') * (~x) ^ ~(mn')) ^ ~q * (~e + ~(f') * (~x) ^ ~(n')) ^ ~(r'), ~x) => (((~x) ^ (~n * FracPart(~q)) * (~c + ~d * (~x) ^ -(~n)) ^ FracPart(~q)) / (~d + ~c * (~x) ^ ~n) ^ FracPart(~q)) * Antiderivative(((~a + ~b * (~x) ^ ~n) ^ ~p * (~d + ~c * (~x) ^ ~n) ^ ~q * (~e + ~f * (~x) ^ ~n) ^ ~r) / (~x) ^ (~n * ~q), ~x) <-- FreeQ([~a, ~b, ~c, ~d, ~e, ~f, ~n, ~p, ~q, ~r], ~x) && (EqQ(~mn, -(~n)) && Not(IntegerQ(~q)))
    @apply_utils Antiderivative((~a + ~(b') * (~x) ^ ~n) ^ ~(p') * (~c + ~(d') * (~x) ^ ~n) ^ ~(q') * (~e1 + ~(f1') * (~x) ^ ~(n2')) ^ ~(r') * (~e2 + ~(f2') * (~x) ^ ~(n2')) ^ ~(r'), ~x) => Antiderivative((~a + ~b * (~x) ^ ~n) ^ ~p * (~c + ~d * (~x) ^ ~n) ^ ~q * (~e1 * ~e2 + ~f1 * ~f2 * (~x) ^ ~n) ^ ~r, ~x) <-- FreeQ([~a, ~b, ~c, ~d, ~e1, ~f1, ~e2, ~f2, ~n, ~p, ~q, ~r], ~x) && (EqQ(~n2, ~n / 2) && (EqQ(~e2 * ~f1 + ~e1 * ~f2, 0) && (IntegerQ(~r) || GtQ(~e1, 0) && GtQ(~e2, 0))))
    @apply_utils Antiderivative((~a + ~(b') * (~x) ^ ~n) ^ ~(p') * (~c + ~(d') * (~x) ^ ~n) ^ ~(q') * (~e1 + ~(f1') * (~x) ^ ~(n2')) ^ ~(r') * (~e2 + ~(f2') * (~x) ^ ~(n2')) ^ ~(r'), ~x) => (((~e1 + ~f1 * (~x) ^ (~n / 2)) ^ FracPart(~r) * (~e2 + ~f2 * (~x) ^ (~n / 2)) ^ FracPart(~r)) / (~e1 * ~e2 + ~f1 * ~f2 * (~x) ^ ~n) ^ FracPart(~r)) * Antiderivative((~a + ~b * (~x) ^ ~n) ^ ~p * (~c + ~d * (~x) ^ ~n) ^ ~q * (~e1 * ~e2 + ~f1 * ~f2 * (~x) ^ ~n) ^ ~r, ~x) <-- FreeQ([~a, ~b, ~c, ~d, ~e1, ~f1, ~e2, ~f2, ~n, ~p, ~q, ~r], ~x) && (EqQ(~n2, ~n / 2) && EqQ(~e2 * ~f1 + ~e1 * ~f2, 0))
end

max_stars_repo_path: src/rules/1 Algebraic functions/1.1 Binomial products/1.1.3 General/1.1.3.5 (a+b x^n)^p (c+d x^n)^q (e+f x^n)^r.jl


text:
#=
Author : Ng Yeesian
Modified : Joshua Ong / Guy Farmer
generates the stochastic deployment model
=#
struct StochasticDeployment_rem <: DeploymentModel
    m::JuMP.Model
    x::Vector{JuMP.VariableRef}
    y::Array{JuMP.VariableRef,3}
    slack::Vector{JuMP.VariableRef}
end
deployment(m::StochasticDeployment_rem) = [round(Int,x) for x in JuMP.value.(m.x)]

function StochasticDeployment_rem(p::DeploymentProblem; nperiods=params.nperiods, tol=params.δ,
    solver=Gurobi.Optimizer(OutputFlag=0), extra_amb = 1)

    old_x = [1, 1, 1, 2, 0, 0, 1, 2, 2, 0, 2, 0, 0, 1, 2, 2, 1, 2, 0, 1, 0, 1, 0, 1, 0, 1, 1, 1, 1, 2, 1, 1, 1, 1, 0, 1, 0, 0, 1, 1, 1, 0, 2, 1]
    nperiods = min(length(p.train), nperiods)
    demand = p.demand[p.train,:]
    I = 1:p.nlocations
    J = 1:p.nregions
    T = 1:nperiods

    m = Model(GLPK.Optimizer) #using GLPK
    JuMP.@variable(m, x[1:p.nlocations] >= 0, Int)
    JuMP.@variable(m, y[1:p.nlocations,1:p.nregions,1:nperiods] >= 0, Int)
    JuMP.@variable(m, z[1:p.nregions,1:nperiods] >= 0, Int)
    JuMP.@variable(m, slack[1:p.nlocations] <= 0, Int) #slack variable

    JuMP.@objective(m, Min, sum(z[j,t] for j=J, t=T) + tol*sum(y[i,j,t] for i=I, j=J, t=T));

    JuMP.@constraint(m, sum(x[i] for i=I) <= p.nambulances - extra_amb);

    #ensure new solution is within slack of the old x.
    for i in I #collect(1:44)
        JuMP.@constraint(m, x[i] - slack[i] == old_x[i]);
    end
    JuMP.@constraint(m, sum(slack[i] for i in I) == extra_amb) #make sure slack is the +/- n ambulances of the old solution

    #Guarantees coverage for all regions
    for j in J
        JuMP.@constraint(m, sum(x[i] for i in filter(i->p.coverage[j,i], I)) >= 1);
    end

    # flow constraints at each station
    for i in I, t in T
        outflow = JuMP.@expression(m, sum(y[i,j,t] for j in filter(j->p.coverage[j,i], J)));
        JuMP.@constraint(m, x[i] >= outflow);
    end

    # shortfall from satisfying demand/calls
    for j in J, t in T
        inflow = JuMP.@expression(m, sum(y[i,j,t] for i in filter(i->p.coverage[j,i], I)));
        JuMP.@constraint(m, z[j,t] >= demand[t,j] - inflow);
    end

    # any location can only house so many hospitals
    for i in I
        JuMP.@constraint(m, x[i] <= 5);
    end

    StochasticDeployment_rem(m, x, y, slack)
end

#status = JuMP.optimize!(model.m)
optimize!(model::StochasticDeployment_rem) = JuMP.optimize!(model.m)
optimize!(model::StochasticDeployment_rem, p::DeploymentProblem) = JuMP.optimize!(model.m)

max_stars_repo_name: michaelhilborn/AmbulanceDeployement

max_issues_repo_name: joshua-ong/AmbulanceDeployment

max_forks_repo_name: joshua-ong/AmbulanceDeployment


text:
n = 1024

while n > 0
    println(n)
    n >>= 1
end

max_stars_repo_name: ethansaxenian/RosettaDecode

max_issues_repo_name: ethansaxenian/RosettaDecode

max_forks_repo_name: ethansaxenian/RosettaDecode


text:
module TimedTests
export @timedtestset, TimedTestSet

import Test
using Test: AbstractTestSet, DefaultTestSet, Broken, Fail, Error, Pass, TestSetException
import Test: record, finish, print_test_errors, print_test_results, print_counts,
            get_testset, get_testset_depth,get_test_counts, get_alignment, filter_errors

macro timedtestset(ex...)
    if length(ex) == 2
        name = ex[1]
        testset = ex[2]
    else
        name = "timed test set"
        testset = ex[1]
    end
    timedtestsetvar = gensym()
    # for some reason, you cannot do @testset TestExtras.TimedTests.TimedTestSet begin ...
    # instead we do: var = TestExtras.TimedTests.TimedTestSet; @testset var begin ...
    return esc(Expr(:block,
    Expr(:(=), timedtestsetvar, :($TimedTestSet)),
    Expr(:macrocall, Symbol("@testset"), __source__, timedtestsetvar, name, testset)))
end

@nospecialize

mutable struct TimedTestSet <: AbstractTestSet
    description::AbstractString
    results::Vector
    n_passed::Int
    anynonpass::Bool
    ti::Float64
    tf::Float64
    TimedTestSet(desc) = new(desc, [], 0, false, time())
end

# For a broken result, simply store the result
Test.record(ts::TimedTestSet, t::Broken) = (push!(ts.results, t); t)
# For a passed result, do not store the result since it uses a lot of memory
record(ts::TimedTestSet, t::Pass) = (ts.n_passed += 1; t)

# For the other result types, immediately print the error message
# but do not terminate. Print a backtrace.
function record(ts::TimedTestSet, t::Union{Fail, Error})
    if Test.TESTSET_PRINT_ENABLE[]
        printstyled(ts.description, ": ", color=:white)
        # don't print for interrupted tests
        if !(t isa Error) || t.test_type !== :test_interrupted
            print(t)
            # don't print the backtrace for Errors because it gets printed in the show
            # method
            if !isa(t, Error)
                Base.show_backtrace(stdout, Test.scrub_backtrace(backtrace()))
            end
            println()
        end
    end
    push!(ts.results, t)
    t, isa(t, Error) || backtrace()
end

# When a TimedTestSet finishes, it records itself to its parent
# testset, if there is one. This allows for recursive printing of
# the results at the end of the tests
record(ts::TimedTestSet, t::AbstractTestSet) = push!(ts.results, t)

@specialize

function print_test_errors(ts::TimedTestSet)
    for t in ts.results
        if isa(t, Error) || isa(t, Fail)
            println("Error in testset $(ts.description):")
            show(t)
            println()
        elseif isa(t, TimedTestSet)
            print_test_errors(t)
        end
    end
end

function print_test_results(ts::TimedTestSet, depth_pad=0)
    # Calculate the overall number for each type so each of
    # the test result types are aligned
    passes, fails, errors, broken, c_passes, c_fails, c_errors, c_broken = get_test_counts(ts)
    total_pass   = passes + c_passes
    total_fail   = fails  + c_fails
    total_error  = errors + c_errors
    total_broken = broken + c_broken
    dig_pass   = total_pass   > 0 ? ndigits(total_pass)   : 0
    dig_fail   = total_fail   > 0 ? ndigits(total_fail)   : 0
    dig_error  = total_error  > 0 ? ndigits(total_error)  : 0
    dig_broken = total_broken > 0 ? ndigits(total_broken) : 0
    total = total_pass + total_fail + total_error + total_broken
    dig_total = total > 0 ? ndigits(total) : 0
    # For each category, take max of digits and header width if there are
    # tests of that type
    pass_width   = dig_pass   > 0 ? max(length("Pass"),   dig_pass)   : 0
    fail_width   = dig_fail   > 0 ? max(length("Fail"),   dig_fail)   : 0
    error_width  = dig_error  > 0 ? max(length("Error"),  dig_error)  : 0
    broken_width = dig_broken > 0 ? max(length("Broken"), dig_broken) : 0
    total_width  = dig_total  > 0 ? max(length("Total"),  dig_total)  : 0
    # Calculate the alignment of the test result counts by
    # recursively walking the tree of test sets
    align = max(get_alignment(ts, 0), length("Test Summary:"))
    # Print the outer test set header once
    pad = total == 0 ? "" : " "
    printstyled(rpad("Test Summary:", align, " "), " |", pad; bold=true, color=:white)
    if pass_width > 0
        printstyled(lpad("Pass", pass_width, " "), "  "; bold=true, color=:green)
    end
    if fail_width > 0
        printstyled(lpad("Fail", fail_width, " "), "  "; bold=true, color=Base.error_color())
    end
    if error_width > 0
        printstyled(lpad("Error", error_width, " "), "  "; bold=true, color=Base.error_color())
    end
    if broken_width > 0
        printstyled(lpad("Broken", broken_width, " "), "  "; bold=true, color=Base.warn_color())
    end
    if total_width > 0
        printstyled(lpad("Total", total_width, " "); bold=true, color=Base.info_color())
    end
    printstyled("  Time: "; bold = true, color = Base.info_color())
    printstyled(string(round(ts.tf-ts.ti; sigdigits=3)); color = Base.info_color())
    println()
    # Recursively print a summary at every level
    print_counts(ts, depth_pad, align, pass_width, fail_width, error_width, broken_width, total_width)
end

# Called at the end of a @testset, behaviour depends on whether
# this is a child of another testset, or the "root" testset
function finish(ts::TimedTestSet)
    # If we are a nested test set, do not print a full summary
    # now - let the parent test set do the printing
    ts.tf = time()
    if get_testset_depth() != 0
        # Attach this test set to the parent test set
        parent_ts = get_testset()
        record(parent_ts, ts)
        return ts
    end
    passes, fails, errors, broken, c_passes, c_fails, c_errors, c_broken = get_test_counts(ts)
    total_pass   = passes + c_passes
    total_fail   = fails  + c_fails
    total_error  = errors + c_errors
    total_broken = broken + c_broken
    total = total_pass + total_fail + total_error + total_broken

    if Test.TESTSET_PRINT_ENABLE[]
        print_test_results(ts)
    end

    # Finally throw an error as we are the outermost test set
    if total != total_pass + total_broken
        # Get all the error/failures and bring them along for the ride
        efs = filter_errors(ts)
        throw(TestSetException(total_pass, total_fail, total_error, total_broken, efs))
    end

    # return the testset so it is returned from the @testset macro
    ts
end

# Recursive function that finds the column that the result counts
# can begin at by taking into account the width of the descriptions
# and the amount of indentation. If a test set had no failures, and
# no failures in child test sets, there is no need to include those
# in calculating the alignment
function get_alignment(ts::TimedTestSet, depth::Int)
    # The minimum width at this depth is
    ts_width = 2*depth + length(ts.description)
    # If all passing, no need to look at children
    !ts.anynonpass && return ts_width
    # Return the maximum of this width and the minimum width
    # for all children (if they exist)
    isempty(ts.results) && return ts_width
    child_widths = map(t->get_alignment(t, depth+1), ts.results)
    return max(ts_width, maximum(child_widths))
end

# Recursive function that fetches backtraces for any and all errors
# or failures the testset and its children encountered
function filter_errors(ts::TimedTestSet)
    efs = []
    for t in ts.results
        if isa(t, Union{TimedTestSet,DefaultTestSet})
            append!(efs, filter_errors(t))
        elseif isa(t, Union{Fail, Error})
            append!(efs, [t])
        end
    end
    efs
end

# Recursive function that counts the number of test results of each
# type directly in the testset, and totals across the child testsets
function get_test_counts(ts::TimedTestSet)
    passes, fails, errors, broken = ts.n_passed, 0, 0, 0
    c_passes, c_fails, c_errors, c_broken = 0, 0, 0, 0
    for t in ts.results
        isa(t, Fail)   && (fails  += 1)
        isa(t, Error)  && (errors += 1)
        isa(t, Broken) && (broken += 1)
        if isa(t, Union{TimedTestSet,DefaultTestSet})
            np, nf, ne, nb, ncp, ncf, nce , ncb = get_test_counts(t)
            c_passes += np + ncp
            c_fails  += nf + ncf
            c_errors += ne + nce
            c_broken += nb + ncb
        end
    end
    ts.anynonpass = (fails + errors + c_fails + c_errors > 0)
    return passes, fails, errors, broken, c_passes, c_fails, c_errors, c_broken
end

# Recursive function that prints out the results at each level of
# the tree of test sets
function print_counts(ts::TimedTestSet, depth, align,
                      pass_width, fail_width, error_width, broken_width, total_width)
    # Count results by each type at this level, and recursively
    # through any child test sets
    passes, fails, errors, broken, c_passes, c_fails, c_errors, c_broken = get_test_counts(ts)
    subtotal = passes + fails + errors + broken + c_passes + c_fails + c_errors + c_broken
    # Print test set header, with an alignment that ensures all
    # the test results appear above each other
    print(rpad(string("  "^depth, ts.description), align, " "), " | ")

    np = passes + c_passes
    if np > 0
        printstyled(lpad(string(np), pass_width, " "), "  ", color=:green)
    elseif pass_width > 0
        # No passes at this level, but some at another level
        print(lpad(" ", pass_width), "  ")
    end

    nf = fails + c_fails
    if nf > 0
        printstyled(lpad(string(nf), fail_width, " "), "  ", color=Base.error_color())
    elseif fail_width > 0
        # No fails at this level, but some at another level
        print(lpad(" ", fail_width), "  ")
    end

    ne = errors + c_errors
    if ne > 0
        printstyled(lpad(string(ne), error_width, " "), "  ", color=Base.error_color())
    elseif error_width > 0
        # No errors at this level, but some at another level
        print(lpad(" ", error_width), "  ")
    end

    nb = broken + c_broken
    if nb > 0
        printstyled(lpad(string(nb), broken_width, " "), "  ", color=Base.warn_color())
    elseif broken_width > 0
        # None broken at this level, but some at another level
        print(lpad(" ", broken_width), "  ")
    end

    if np == 0 && nf == 0 && ne == 0 && nb == 0
        printstyled("No tests", color=Base.info_color())
    else
        printstyled(lpad(string(subtotal), total_width, " "), color=Base.info_color())
    end
    println()

    # Only print results at lower levels if we had failures
    if np + nb != subtotal
        for t in ts.results
            if isa(t, Union{TimedTestSet,DefaultTestSet})
                print_counts(t, depth + 1, align,
                    pass_width, fail_width, error_width, broken_width, total_width)
            end
        end
    end
end

end


text:
# --------------------------------------------------------------------------- #
# collect the un-hidden filenames available in a given directory

function getfname(target)
    # target : string := chemin + nom du repertoire ou se trouve les instances

    # positionne le currentdirectory dans le repertoire cible
    cd(joinpath(homedir(),target))

    # retourne le repertoire courant
    println("pwd = ", pwd())

    # recupere tous les fichiers se trouvant dans le repertoire data
    allfiles = readdir()

    # vecteur booleen qui marque les noms de fichiers valides
    flag = trues(size(allfiles))

    k=1
    for f in allfiles
        # traite chaque fichier du repertoire
        if f[1] != '.'
            # pas un fichier cache => conserver
            println("fname = ", f)
        else
            # fichier cache => supprimer
            flag[k] = false
        end
        k = k+1
    end

    # extrait les noms valides et retourne le vecteur correspondant
    finstances = allfiles[flag]
    return finstances
end