text:
"""
    MoistThermodynamics

Moist thermodynamic functions, e.g., for air pressure (atmosphere equation
of state), latent heats of phase transitions, saturation vapor pressures, and
saturation specific humidities
"""
module MoistThermodynamics

using DocStringExtensions

using ..RootSolvers
using ...PlanetParameters

# Atmospheric equation of state
export air_pressure, air_temperature, air_density, specific_volume, soundspeed_air

# Energies
export total_energy, internal_energy, internal_energy_sat

# Specific heats of moist air
export cp_m, cv_m, gas_constant_air, moist_gas_constants

# Latent heats
export latent_heat_vapor, latent_heat_sublim, latent_heat_fusion

# Saturation vapor pressures and specific humidities over liquid and ice
export Liquid, Ice
export saturation_vapor_pressure, saturation_shum_generic, saturation_shum
export saturation_excess

# Functions used in thermodynamic equilibrium among phases (liquid and ice
# determined diagnostically from total water specific humidity)
export liquid_fraction_equil, liquid_fraction_nonequil, saturation_adjustment, PhasePartition_equil

# Auxiliary functions, e.g., for diagnostic purposes
export air_temperature_from_liquid_ice_pottemp, dry_pottemp, virtual_pottemp, exner
export liquid_ice_pottemp, liquid_ice_pottemp_sat

include("states.jl")

"""
    gas_constant_air([q::PhasePartition])

The specific gas constant of moist air given
 - `q` [`PhasePartition`](@ref). Without this argument the results are that of dry air.
"""
gas_constant_air(q::PhasePartition{DT}) where {DT} =
  DT(R_d) * ( 1 +  (DT(molmass_ratio) - 1)*q.tot - DT(molmass_ratio)*(q.liq + q.ice) )
gas_constant_air() = gas_constant_air(PhasePartition(0.0))

"""
    gas_constant_air(ts::ThermodynamicState)

The specific gas constant of moist air given
a thermodynamic state `ts`.
"""
gas_constant_air(ts::ThermodynamicState) =
  gas_constant_air(PhasePartition(ts))

"""
    air_pressure(T, ρ[, q::PhasePartition])

The air pressure from the equation of state
(ideal gas law) where:

 - `T` air temperature
 - `ρ` (moist-)air density
and, optionally,
 - `q` [`PhasePartition`](@ref). Without this argument the results are that of dry air.
"""
air_pressure(T::DT, ρ::DT, q::PhasePartition=PhasePartition(zero(DT))) where {DT<:Real} =
  gas_constant_air(q) * ρ * T

"""
    air_pressure(ts::ThermodynamicState)

The air pressure from the equation of state
(ideal gas law), given a thermodynamic state `ts`.
"""
air_pressure(ts::ThermodynamicState) =
  air_pressure(air_temperature(ts), air_density(ts), PhasePartition(ts))

"""
    air_density(T, p[, q::PhasePartition])

The (moist-)air density from the equation of state
(ideal gas law) where

 - `T` air temperature
 - `p` pressure
and, optionally,
 - `q` [`PhasePartition`](@ref). Without this argument the results are that of dry air.
"""
air_density(T::DT, p::DT, q::PhasePartition=PhasePartition(zero(DT))) where {DT<:Real} =
  p / (gas_constant_air(q) * T)

"""
    air_density(ts::ThermodynamicState)

The (moist-)air density from the equation of state
(ideal gas law), given a thermodynamic state `ts`.
"""
air_density(ts::ThermodynamicState) = ts.density

"""
    specific_volume(T, p[, q::PhasePartition])

The (moist-)air specific volume from the equation of
state (ideal gas law) where

 - `T` air temperature
 - `p` pressure
and, optionally,
 - `q` [`PhasePartition`](@ref). Without this argument the results are that of dry air.
"""
specific_volume(T::DT, p::DT, q::PhasePartition=PhasePartition(zero(DT))) where {DT<:Real} =
  (gas_constant_air(q) * T) / p

"""
    specific_volume(ts::ThermodynamicState)

The (moist-)air specific volume from the equation of
state (ideal gas law), given a thermodynamic state `ts`.
"""
specific_volume(ts::ThermodynamicState) = specific_volume(air_temperature(ts), air_pressure(ts), PhasePartition(ts))

"""
    cp_m([q::PhasePartition])

The isobaric specific heat capacity of moist
air where, optionally,
 - `q` [`PhasePartition`](@ref). Without this argument the results are that of dry air.
"""
cp_m(q::PhasePartition{DT}) where {DT} =
  DT(cp_d) + (DT(cp_v) - DT(cp_d))*q.tot + (DT(cp_l) - DT(cp_v))*q.liq + (DT(cp_i) - DT(cp_v))*q.ice
cp_m() = cp_m(PhasePartition(0.0))

"""
    cp_m(ts::ThermodynamicState)

The isobaric specific heat capacity of moist
air, given a thermodynamic state `ts`.
"""
cp_m(ts::ThermodynamicState) = cp_m(PhasePartition(ts))

"""
    cv_m([q::PhasePartition])

The isochoric specific heat capacity of moist
air where optionally,
 - `q` [`PhasePartition`](@ref). Without this argument the results are that of dry air.
"""
cv_m(q::PhasePartition{DT}) where {DT} =
  DT(cv_d) + (DT(cv_v) - DT(cv_d))*q.tot + (DT(cv_l) - DT(cv_v))*q.liq + (DT(cv_i) - DT(cv_v))*q.ice
cv_m() = cv_m(PhasePartition(0.0))

"""
    cv_m(ts::ThermodynamicState)

The isochoric specific heat capacity of moist
air given a thermodynamic state `ts`.
"""
cv_m(ts::ThermodynamicState) = cv_m(PhasePartition(ts))

"""
    (R_m, cp_m, cv_m, γ_m) = moist_gas_constants([q::PhasePartition])

Wrapper to compute all gas constants at once, where optionally,
 - `q` [`PhasePartition`](@ref). Without this argument the results are that of dry air.

The function returns a tuple of
 - `R_m` [`gas_constant_air`](@ref)
 - `cp_m` [`cp_m`](@ref)
 - `cv_m` [`cv_m`](@ref)
 - `γ_m = cp_m/cv_m`

Without the specific humidity arguments, the results
are that of dry air.
"""
function moist_gas_constants(q::PhasePartition=PhasePartition(0.0))
    R_gas  = gas_constant_air(q)
    cp = cp_m(q)
    cv = cv_m(q)
    γ = cp/cv

    return (R_gas, cp, cv, γ)
end

"""
    (R_m, cp_m, cv_m, γ_m) = moist_gas_constants(ts::ThermodynamicState)

Wrapper to compute all gas constants at once, given a thermodynamic state `ts`.

The function returns a tuple of
 - `R_m` [`gas_constant_air`](@ref)
 - `cp_m` [`cp_m`](@ref)
 - `cv_m` [`cv_m`](@ref)
 - `γ_m = cp_m/cv_m`

"""
moist_gas_constants(ts::ThermodynamicState) = moist_gas_constants(PhasePartition(ts))

"""
    air_temperature(e_int, q::PhasePartition)

The air temperature, where

 - `e_int` internal energy per unit mass
and, optionally,
 - `q` [`PhasePartition`](@ref). Without this argument the results are that of dry air.
"""
function air_temperature(e_int::DT, q::PhasePartition=PhasePartition(zero(DT))) where {DT}
  T_0 +
    (e_int - (q.tot - q.liq) * DT(e_int_v0) + q.ice * (DT(e_int_v0) + DT(e_int_i0))) /
    cv_m(q)
end

"""
    air_temperature(ts::ThermodynamicState)

The air temperature, given a thermodynamic state `ts`.
"""
air_temperature(ts::PhaseEquil) = ts.T
air_temperature(ts::PhaseNonEquil) = air_temperature(ts.e_int, PhasePartition(ts))

"""
    internal_energy(T[, q::PhasePartition])

The internal energy per unit mass, given a thermodynamic state `ts` or

 - `T` temperature
and, optionally,
 - `q` [`PhasePartition`](@ref). Without this argument the results are that of dry air.
"""
internal_energy(T::DT, q::PhasePartition=PhasePartition(zero(DT))) where {DT} =
  cv_m(q) * (T - DT(T_0)) +
  (q.tot - q.liq) * DT(e_int_v0) - q.ice * (DT(e_int_v0) + DT(e_int_i0))

"""
    internal_energy(ts::ThermodynamicState)

The internal energy per unit mass, given a thermodynamic state `ts`.
"""
internal_energy(ts::PhaseEquil) = ts.e_int
internal_energy(ts::PhaseNonEquil) = ts.e_int

"""
    internal_energy_sat(T, ρ, q_tot)

The internal energy per unit mass in thermodynamic equilibrium at saturation where

 - `T` temperature
 - `ρ` (moist-)air density
 - `q_tot` total specific humidity
"""
internal_energy_sat(T::DT, ρ::DT, q_tot::DT) where DT =
  internal_energy(T, PhasePartition_equil(T, ρ, q_tot))

"""
    internal_energy_sat(ts::ThermodynamicState)

The internal energy per unit mass in
thermodynamic equilibrium at saturation,
given a thermodynamic state `ts`.
"""
internal_energy_sat(ts::PhaseEquil) =
  internal_energy_sat(air_temperature(ts), air_density(ts), ts.q_tot)
internal_energy_sat(ts::PhaseNonEquil) =
  internal_energy_sat(air_temperature(ts), air_density(ts), ts.q.tot)

"""
    total_energy(e_kin, e_pot, T[, q::PhasePartition])

The total energy per unit mass, given

 - `e_kin` kinetic energy per unit mass
 - `e_pot` potential energy per unit mass
 - `T` temperature
and, optionally,
 - `q` [`PhasePartition`](@ref). Without this argument the results are that of dry air.

"""
total_energy(e_kin::DT, e_pot::DT, T::DT, q::PhasePartition=PhasePartition(zero(DT))) where {DT} =
  e_kin + e_pot + internal_energy(T, q)

"""
    total_energy(e_kin, e_pot, ts::ThermodynamicState)

The total energy per unit mass
given a thermodynamic state `ts`.
"""
total_energy(e_kin, e_pot, ts::ThermodynamicState) = internal_energy(ts) + e_kin + e_pot

"""
    soundspeed_air(T[, q::PhasePartition])

The speed of sound in air, where
 - `T` temperature
and, optionally,
 - `q` [`PhasePartition`](@ref). Without this argument the results are that of dry air.

Without the specific humidity arguments, the results
are that of dry air.
"""
function soundspeed_air(T::DT, q::PhasePartition=PhasePartition(zero(DT))) where {DT<:Real}
  γ   = cp_m(q) / cv_m(q)
  R_m = gas_constant_air(q)
  return sqrt(γ*R_m*T)
end

"""
    soundspeed_air(ts::ThermodynamicState)

The speed of sound in air given a thermodynamic state `ts`.
"""
soundspeed_air(ts::ThermodynamicState) =
  soundspeed_air(air_temperature(ts), PhasePartition(ts))

"""
    latent_heat_vapor(T::Real)

The specific latent heat of vaporization where
 - `T` temperature
"""
latent_heat_vapor(T::DT) where {DT<:Real} =
  latent_heat_generic(T, DT(LH_v0), DT(cp_v) - DT(cp_l))

"""
    latent_heat_vapor(ts::ThermodynamicState)

The specific latent heat of vaporization
given a thermodynamic state `ts`.
"""
latent_heat_vapor(ts::ThermodynamicState{DT}) where {DT} =
    latent_heat_vapor(air_temperature(ts))

"""
    latent_heat_sublim(T::Real)

The specific latent heat of sublimation where
 - `T` temperature
"""
latent_heat_sublim(T::DT) where {DT<:Real} =
  latent_heat_generic(T, DT(LH_s0), DT(cp_v) - DT(cp_i))

"""
    latent_heat_sublim(ts::ThermodynamicState)

The specific latent heat of sublimation
given a thermodynamic state `ts`.
"""
latent_heat_sublim(ts::ThermodynamicState{DT}) where {DT} =
  latent_heat_sublim(air_temperature(ts))

"""
    latent_heat_fusion(T)

The specific latent heat of fusion where
 - `T` temperature
"""
latent_heat_fusion(T::DT) where {DT<:Real} =
  latent_heat_generic(T, DT(LH_f0), DT(cp_l) - DT(cp_i))

"""
    latent_heat_fusion(ts::ThermodynamicState)

The specific latent heat of fusion
given a thermodynamic state `ts`.
"""
latent_heat_fusion(ts::ThermodynamicState{DT}) where {DT} =
  latent_heat_fusion(air_temperature(ts))

"""
    latent_heat_generic(T::Real, LH_0::Real, Δcp::Real)

Return the specific latent heat of a generic phase transition between
two phases using Kirchhoff's relation.

The latent heat computation assumes constant isobaric specifc heat capacities
of the two phases. `T` is the temperature, `LH_0` is the latent heat of the
phase transition at `T_0`, and `Δcp` is the difference between the isobaric
specific heat capacities (heat capacity in the higher-temperature phase minus
that in the lower-temperature phase).
"""
latent_heat_generic(T::Real, LH_0::Real, Δcp::Real) =
  LH_0 + Δcp * (T - T_0)

"""
    Phase

A phase condensate, to dispatch over
[`saturation_vapor_pressure`](@ref) and
[`saturation_shum_generic`](@ref).
"""
abstract type Phase end

"""
    Liquid <: Phase

A liquid phase, to dispatch over
[`saturation_vapor_pressure`](@ref) and
[`saturation_shum_generic`](@ref).
"""
struct Liquid <: Phase end

"""
    Ice <: Phase

An ice phase, to dispatch over
[`saturation_vapor_pressure`](@ref) and
[`saturation_shum_generic`](@ref).
"""
struct Ice <: Phase end

"""
    saturation_vapor_pressure(T, Liquid())

Return the saturation vapor pressure over a plane liquid surface at
temperature `T`.

    saturation_vapor_pressure(T, Ice())

Return the saturation vapor pressure over a plane ice surface at
temperature `T`.

    saturation_vapor_pressure(T, LH_0, Δcp)

Compute the saturation vapor pressure over a plane surface by integration
of the Clausius-Clepeyron relation.

The Clausius-Clapeyron relation

    dlog(p_v_sat)/dT = [LH_0 + Δcp * (T-T_0)]/(R_v*T^2)

is integrated from the triple point temperature `T_triple`, using
Kirchhoff's relation

    L = LH_0 + Δcp * (T - T_0)

for the specific latent heat `L` with constant isobaric specific
heats of the phases. The linear dependence of the specific latent heat
on temperature `T` allows analytic integration of the Clausius-Clapeyron
relation to obtain the saturation vapor pressure `p_v_sat` as a function of
the triple point pressure `press_triple`.

"""
saturation_vapor_pressure(T::DT, ::Liquid) where DT = saturation_vapor_pressure(T, DT(LH_v0), DT(cp_v) - DT(cp_l))
function saturation_vapor_pressure(ts::ThermodynamicState, ::Liquid)

    return saturation_vapor_pressure(air_temperature(ts), LH_v0, cp_v - cp_l)

end
saturation_vapor_pressure(T::DT, ::Ice) where DT = saturation_vapor_pressure(T, DT(LH_s0), DT(cp_v) - DT(cp_i))
saturation_vapor_pressure(ts::ThermodynamicState, ::Ice) = saturation_vapor_pressure(air_temperature(ts), LH_s0, cp_v - cp_i)

function saturation_vapor_pressure(T::DT, LH_0::DT, Δcp::DT) where DT

    return DT(press_triple) * (T/DT(T_triple))^(Δcp/DT(R_v)) *
        exp( (DT(LH_0) - Δcp*DT(T_0))/DT(R_v) * (1 / DT(T_triple) - 1 / T) )

end

"""
    saturation_shum_generic(T, ρ[; phase=Liquid()])

Compute the saturation specific humidity over a plane surface of condensate, given

 - `T` temperature
 - `ρ` (moist-)air density
and, optionally,
 - `Liquid()` indicating condensate is liquid
 - `Ice()` indicating condensate is ice
"""
function saturation_shum_generic(T::DT, ρ::DT; phase::Phase=Liquid()) where DT
    p_v_sat = saturation_vapor_pressure(T, phase)
    return saturation_shum_from_pressure(T, ρ, p_v_sat)
end

"""
    saturation_shum(T, ρ[, q::PhasePartition])

Compute the saturation specific humidity, given

 - `T` temperature
 - `ρ` (moist-)air density
and, optionally,
 - `q` [`PhasePartition`](@ref)

If the phase partition `q` is given, the saturation specific humidity is that over a
mixture of liquid and ice, computed in a thermodynamically consistent way from the
weighted sum of the latent heats of the respective phase transitions (Pressel et al.,
JAMES, 2015).  That is, the saturation vapor pressure and from it the saturation specific
humidity are computed from a weighted mean of the latent heats of vaporization and
sublimation, with the weights given by the fractions of condensate `q.liq/(q.liq + q.ice)`
and `q.ice/(q.liq + q.ice)` that are liquid and ice, respectively.

If the condensate specific humidities `q.liq` and `q.ice` are not given or are both
zero, the saturation specific humidity is that over a mixture of liquid and ice,
with the fraction of liquid given by temperature dependent `liquid_fraction_equil(T)`
and the fraction of ice by the complement `1 - liquid_fraction_equil(T)`.
"""
function saturation_shum(T::DT, ρ::DT, q::PhasePartition=PhasePartition(zero(DT))) where DT

    # get phase partitioning
    _liquid_frac = liquid_fraction_equil(T, q)
    _ice_frac    = 1 - _liquid_frac

    # effective latent heat at T_0 and effective difference in isobaric specific
    # heats of the mixture
    LH_0    = _liquid_frac * DT(LH_v0) + _ice_frac * DT(LH_s0)
    Δcp     = _liquid_frac * (DT(cp_v) - DT(cp_l)) + _ice_frac * (DT(cp_v) - DT(cp_i))

    # saturation vapor pressure over possible mixture of liquid and ice
    p_v_sat = saturation_vapor_pressure(T, DT(LH_0), Δcp)

    return saturation_shum_from_pressure(T, ρ, p_v_sat)

end

"""
    saturation_shum(ts::ThermodynamicState)

Compute the saturation specific humidity, given a thermodynamic state `ts`.
"""
saturation_shum(ts::ThermodynamicState) =
  saturation_shum(air_temperature(ts), air_density(ts), PhasePartition(ts))

"""
    saturation_shum_from_pressure(T, ρ, p_v_sat)

Compute the saturation specific humidity, given

 - `T` ambient air temperature,
 - `ρ` density
 - `p_v_sat` saturation vapor pressure
"""
saturation_shum_from_pressure(T::DT, ρ::DT, p_v_sat::DT) where {DT<:Real} =
  min(1, p_v_sat / (ρ * DT(R_v) * T))

"""
    saturation_excess(T, ρ, q::PhasePartition)

The saturation excess in equilibrium where

 - `T` temperature
 - `ρ` (moist-)air density
 - `q` [`PhasePartition`](@ref)

The saturation excess is the difference between the total specific humidity `q.tot`
and the saturation specific humidity in equilibrium, and it is defined to be
nonzero only if this difference is positive.
"""
saturation_excess(T::DT, ρ::DT, q::PhasePartition) where {DT} =
  max(0, q.tot - saturation_shum(T, ρ, q))

"""
    saturation_excess(ts::ThermodynamicState)

Compute the saturation excess in equilibrium,
given a thermodynamic state `ts`.
"""
saturation_excess(ts::ThermodynamicState) =
  saturation_excess(air_temperature(ts), air_density(ts), PhasePartition(ts))

"""
    liquid_fraction_equil(T[, q::PhasePartition])

The fraction of condensate, assuming phase equilibrium, that is liquid where

 - `T` temperature
 - `q` [`PhasePartition`](@ref)

If `q.liq` or `q.ice` are nonzero, the liquid fraction is computed from
them.

Otherwise, the fraction of liquid is a function that is 1 above `T_freeze` and goes to
zero below `T_freeze`.
"""
function liquid_fraction_equil(T::DT, q::PhasePartition=PhasePartition(zero(DT))) where {DT<:Real}
  q_c = q.liq + q.ice     # condensate specific humidity
  if q_c > 0
    return q.liq / q_c
  else
    # For now: Heaviside function for partitioning into liquid and ice: all liquid
    # for T > T_freeze; all ice for T <= T_freeze
    return DT(T > DT(T_freeze))
  end
end

"""
    liquid_fraction_equil(ts::ThermodynamicState)

The fraction of condensate that is liquid given a thermodynamic state `ts`.
"""
liquid_fraction_equil(ts::ThermodynamicState) =
  liquid_fraction_equil(air_temperature(ts), PhasePartition(ts))

"""
    liquid_fraction_nonequil(T[, q::PhasePartition])

The fraction of condensate, assuming phase non-equilibrium, that is liquid where

 - `T` temperature
 - `q` [`PhasePartition`](@ref)

If `q.liq` or `q.ice` are nonzero, the liquid fraction is computed from
them.

Otherwise, the fraction of liquid is a function that is 1 above `T_freeze` and goes to
zero below `T_freeze`.

!!! todo
    Currently [`liquid_fraction_nonequil`](@ref) calls [`liquid_fraction_equil`](@ref),
    but we should implement a more general function here.
"""
liquid_fraction_nonequil(T::DT, q::PhasePartition) where {DT<:Real} =
  liquid_fraction_equil(T, q)

"""
    liquid_fraction_nonequil(ts::ThermodynamicState)

The fraction of condensate that is liquid given a thermodynamic state `ts`.
"""
liquid_fraction_nonequil(ts::ThermodynamicState) =
  liquid_fraction_nonequil(air_temperature(ts), PhasePartition(ts))

"""
    PhasePartition_equil(T, ρ, q_tot)

Partition the phases in equilibrium, returning a [`PhasePartition`](@ref) object using the
[`liquid_fraction_equil`](@ref) function where

 - `T` temperature
 - `ρ` (moist-)air density
 - `q_tot` total specific humidity

The residual `q.tot - q.liq - q.ice` is the vapor specific humidity.
"""
function PhasePartition_equil(T::DT, ρ::DT, q_tot::DT) where {DT}
    _liquid_frac = liquid_fraction_equil(T)   # fraction of condensate that is liquid
    q_c   = saturation_excess(T, ρ, PhasePartition(q_tot))   # condensate specific humidity
    q_liq = _liquid_frac * q_c  # liquid specific humidity
    q_ice = (1 - _liquid_frac) * q_c # ice specific humidity

    return PhasePartition(q_tot, q_liq, q_ice)
end

PhasePartition(ts::PhaseEquil) = PhasePartition_equil(air_temperature(ts), air_density(ts), ts.q_tot)
PhasePartition(ts::PhaseNonEquil) = ts.q

"""
    saturation_adjustment(e_int, ρ, q_tot)

Compute the temperature that is consistent with

 - `e_int` internal energy
 - `ρ` (moist-)air density
 - `q_tot` total specific humidity

See also [`saturation_adjustment_q_t_θ_l`](@ref).
"""
function saturation_adjustment(e_int::DT, ρ::DT, q_tot::DT) where DT
  T_1 = max(DT(T_min), air_temperature(e_int, PhasePartition(q_tot))) # Assume all vapor
  q_v_sat = saturation_shum(T_1, ρ)
  if q_tot <= q_v_sat # If not saturated return T_1
    return T_1
  else # If saturated, iterate
    # FIXME here: need to revisit bounds for saturation adjustment to guarantee bracketing of zero.
    T_2 = air_temperature(e_int, PhasePartition(q_tot, DT(0), q_tot)) # Assume all ice
    T, converged = find_zero(
      T -> internal_energy_sat(T, ρ, q_tot) - e_int,
      T_1, T_2, SecantMethod(); xatol=DT(1e-3), maxiters=10)
    return T
  end
end

"""
    saturation_adjustment_q_tot_θ_liq_ice(θ_liq_ice, q_tot, ρ, p)

Compute the temperature that is consistent with

 - `θ_liq_ice` liquid-ice potential temperature
 - `q_tot` total specific humidity
 - `ρ` density
 - `p` pressure

See also [`saturation_adjustment`](@ref).
"""
function saturation_adjustment_q_tot_θ_liq_ice(θ_liq_ice::DT, q_tot::DT, ρ::DT, p::DT) where DT
  T_1 = air_temperature_from_liquid_ice_pottemp(θ_liq_ice, p) # Assume all vapor
  q_v_sat = saturation_shum(T_1, ρ)
  if q_tot <= q_v_sat # If not saturated
    return T_1
  else  # If saturated, iterate
    T_2 = air_temperature_from_liquid_ice_pottemp(θ_liq_ice, p, PhasePartition(q_tot, DT(0), q_tot)) # Assume all ice
    T, converged = find_zero(
      T -> θ_liq_ice - liquid_ice_pottemp_sat(T, p, PhasePartition_equil(T, ρ, q_tot)),
      T_1, T_2, SecantMethod(); xatol=DT(1e-3), maxiters=10)
    return T
  end
end

"""
    liquid_ice_pottemp(T, p, q::PhasePartition)

The liquid-ice potential temperature where
 - `T` temperature
 - `p` pressure
and, optionally,
 - `q` [`PhasePartition`](@ref). Without this argument the results are that of dry air.
"""
function liquid_ice_pottemp(T::DT, p::DT, q::PhasePartition=PhasePartition(zero(DT))) where DT
    # isobaric specific heat of moist air
    _cp_m   = cp_m(q)

    # liquid-ice potential temperature, approximating latent heats
    # of phase transitions as constants
    return dry_pottemp(T, p, q) * (1 - (DT(LH_v0)*q.liq + DT(LH_s0)*q.ice)/(_cp_m*T))
end

"""
    liquid_ice_pottemp(ts::ThermodynamicState)

The liquid-ice potential temperature,
given a thermodynamic state `ts`.
"""
liquid_ice_pottemp(ts::ThermodynamicState) =
  liquid_ice_pottemp(air_temperature(ts), air_pressure(ts), PhasePartition(ts))

"""
    dry_pottemp(T, p[, q::PhasePartition])

The dry potential temperature where

 - `T` temperature
 - `p` pressure
and, optionally,
 - `q` [`PhasePartition`](@ref). Without this argument the results are that of dry air.
 """
dry_pottemp(T::DT, p::DT, q::PhasePartition=PhasePartition(zero(DT))) where {DT<:Real} =
  T / exner(p, q)

"""
    dry_pottemp(ts::ThermodynamicState)

The dry potential temperature, given a thermodynamic state `ts`.
"""
dry_pottemp(ts::ThermodynamicState) =
  dry_pottemp(air_temperature(ts), air_pressure(ts), PhasePartition(ts))

"""
    air_temperature_from_liquid_ice_pottemp(θ_liq_ice, p[, q::PhasePartition])

The air temperature, where

 - `θ_liq_ice` liquid-ice potential temperature
 - `p` pressure
and, optionally,
 - `q` [`PhasePartition`](@ref). Without this argument the results are that of dry air.
"""
air_temperature_from_liquid_ice_pottemp(θ_liq_ice::DT, p::DT, q::PhasePartition=PhasePartition(zero(DT))) where {DT} =
  θ_liq_ice*exner(p, q) + (DT(LH_v0)*q.liq + DT(LH_s0)*q.ice) / cp_m(q)

"""
    virtual_pottemp(T, p[, q::PhasePartition])

The virtual temperature where

 - `T` temperature
 - `p` pressure
and, optionally,
 - `q` [`PhasePartition`](@ref). Without this argument the results are that of dry air.
"""
virtual_pottemp(T::DT, p::DT, q::PhasePartition=PhasePartition(zero(DT))) where {DT} =
  gas_constant_air(q) / R_d * dry_pottemp(T, p, q)

"""
    virtual_pottemp(ts::ThermodynamicState)

The virtual potential temperature,
given a thermodynamic state `ts`.
"""
virtual_pottemp(ts::ThermodynamicState) =
  virtual_pottemp(air_temperature(ts), air_pressure(ts), PhasePartition(ts))

"""
    liquid_ice_pottemp_sat(T, p[, q::PhasePartition])

The saturated liquid ice potential temperature where

 - `T` temperature
 - `p` pressure
and, optionally,
 - `q` [`PhasePartition`](@ref). Without this argument the results are that of dry air.
"""
function liquid_ice_pottemp_sat(T::DT, p::DT, q::PhasePartition=PhasePartition(zero(DT))) where {DT}
    ρ = air_density(T, p, q)
    q_v_sat = saturation_shum(T, ρ, q)
    return liquid_ice_pottemp(T, p, PhasePartition(q_v_sat))
end

"""
    liquid_ice_pottemp_sat(ts::ThermodynamicState)

The liquid potential temperature given a thermodynamic state `ts`.
"""
liquid_ice_pottemp_sat(ts::ThermodynamicState) =
  liquid_ice_pottemp_sat(air_temperature(ts), air_pressure(ts), PhasePartition(ts))

"""
    exner(p[, q::PhasePartition])

Compute the Exner function where
 - `p` pressure
and, optionally,
 - `q` [`PhasePartition`](@ref). Without this argument the results are that of dry air.
"""
function exner(p::DT, q::PhasePartition=PhasePartition(zero(DT))) where {DT}
    # gas constant and isobaric specific heat of moist air
    _R_m    = gas_constant_air(q)
    _cp_m   = cp_m(q)

    return (p/DT(MSLP))^(_R_m/_cp_m)
end

"""
    exner(ts::ThermodynamicState)

Compute the Exner function, given a thermodynamic state `ts`.
"""
exner(ts::ThermodynamicState) =
  exner(air_pressure(ts), PhasePartition(ts))

end #module MoistThermodynamics.jl

max_stars_repo_name: christophernhill-dev/CLIMA

max_issues_repo_name: christophernhill-dev/CLIMA

max_forks_repo_name: christophernhill-dev/CLIMA


text:
using ClimateTools
using AxisArrays
using Shapefile
using NetCDF
using NCDatasets
using Dates
using Statistics
using Random
using ArgCheck
using DataFrames
using Test
using HDF5
using Extremes

include("interface_test.jl")
include("utils_test.jl")
include("extract_test.jl")
include("biascorrect_test.jl")
include("indices_test.jl")
# include("mapping_test.jl")
include("functions_test.jl")
include("spatial_test.jl")

ext: jl


text:
# stolen from python implementation
# https://technicallyagarwal.wordpress.com/2018/02/06/markov-chain-algorithm/

mutable struct Table
    prefixVec::Vector{Tuple{String, String}}
    suffixVec::Vector{String}
    Table() = new(Vector{Tuple{String,String}}(), Vector{String}())
end

function Table(s::String)
    splits = split(s)
    buffer = Table()
    for i in eachindex(splits[1:end-2])
        prefixes::Tuple{String, String} = tuple(String(splits[i]),
                                                String(splits[i+1]))
        suffix::String = splits[i+2]
        if !(prefixes in buffer.prefixVec)
            push!(buffer.prefixVec, prefixes)
            push!(buffer.suffixVec, suffix)
        end
    end
    return buffer
end

function lookup(prefix::Tuple{String, String}, table::Table)
   index::Union{Nothing,Int64} = findfirst(x-> prefix == x, table.prefixVec)
    if index != nothing
        return table.suffixVec[index]
    end
    return ""
end

function generate_markov_text(table::Table)
    buffer = ""
    n::Int64 = rand(20:100)
    items = shuffle(table.prefixVec), shuffle(table.suffixVec)
    currentPrefixes, currentSuffixes = items[1][1], items[2][1]
    buffer = @sprintf("%s %s ", currentPrefixes[1], currentPrefixes[2])

    for i in 1:n
        if isempty(currentSuffixes)
            currentPrefixes = rand(items[1])
            currentSuffixes = rand(items[2])
        else
            word = rand(items[2])
            buffer *= @sprintf("%s ", word)
            currentPrefixes = (currentPrefixes[2], word)
            currentSuffixes = lookup(currentPrefixes, table)
        end
    end

    return buffer
end

function generate_markov_from_string(s::String)
    temp = Table(s)
    gen_text = generate_markov_text(temp)
    return gen_text
end

max_stars_repo_name: Iaggelis/bald-bot

max_issues_repo_name: Iaggelis/bald-bot

max_forks_repo_name: Iaggelis/bald-bot


text:
# ------------------------------------------------------------------------------------------
# # Introduction to Recurrent Neural Networks
# (c) Deniz Yuret, 2019
# * Objectives: learn about RNNs, the RNN layer, compare with MLP on a tagging task.
# * Prerequisites: [MLP models](40.mlp.ipynb)
# * New functions:
# [RNN](http://denizyuret.github.io/Knet.jl/latest/reference/#Knet.RNN),
# [adam](http://denizyuret.github.io/Knet.jl/latest/reference/#Knet.adam)
#
# ![image](https://github.com/denizyuret/Knet.jl/blob/master/docs/src/images/RNN-
# unrolled.png?raw=true)([image
# source](http://colah.github.io/posts/2015-08-Understanding-LSTMs))
#
# In this notebook we will see how to implement a recurrent neural network (RNN) in Knet. In
# RNNs, connections between units form a directed cycle, which allows them to keep a
# persistent state over time. This gives them the ability to process sequences of arbitrary
# length one element at a time, while keeping track of what happened at previous elements.
# One can view the current state of the RNN as a representation for the sequence processed
# so far.
#
# We will build a part-of-speech tagger using a large annotated corpus of English. We will
# represent words with numeric vectors appropriate as inputs to a neural network. These word
# vectors will be initialized randomly and learned during training just like other model
# parameters. We will compare three network architectures: (1) an MLP which tags each word
# independently of its neighbors, (2) a simple RNN that can represent the neighboring words
# to the left, (3) a bidirectional RNN that can represent both left and right contexts. As
# can be expected 1 < 2 < 3 in performance. More surprisingly, the three models are very
# similar to each other: we will see their model diagrams are identical except for the
# horizontal connections that carry information across the sequence.
# ------------------------------------------------------------------------------------------

# Setup display width, load packages, import symbols
ENV["COLUMNS"] = 72
using Pkg; for p in ("Knet","Plots"); haskey(Pkg.installed(),p) || Pkg.add(p); end
using Random: shuffle!
using Base.Iterators: flatten
using Knet: Knet, AutoGrad, param, param0, mat, RNN, relu, Data, adam, progress, nll, zeroone

# ------------------------------------------------------------------------------------------
# ## The Brown Corpus
# To introduce recurrent neural networks (RNNs) we will train a part-of-speech tagger using
# the [Brown Corpus](https://en.wikipedia.org/wiki/Brown_Corpus). We will train three
# models: a MLP, a unidirectional RNN, a bidirectional RNN and observe significant
# performance differences.
# ------------------------------------------------------------------------------------------

include(Knet.dir("data/nltk.jl"))
(data,words,tags) = brown()
println("The Brown Corpus has $(length(data)) sentences, $(sum(length(p[1]) for p in data)) tokens, with a word vocabulary of $(length(words)) and a tag vocabulary of $(length(tags)).")

# ------------------------------------------------------------------------------------------
# `data` is an array of `(w,t)` pairs each representing a sentence, where `w` is a sequence
# of word ids, and `t` is a sequence of tag ids. `words` and `tags` contain the strings for
# the ids.
# ------------------------------------------------------------------------------------------

println.(summary.((data,words,tags)));

# ------------------------------------------------------------------------------------------
# Here is what the first sentence looks like with ids and with strings:
# ------------------------------------------------------------------------------------------

(w,t) = first(data)
display(permutedims(Int[w t]))
display(permutedims([words[w] tags[t]]))

# ------------------------------------------------------------------------------------------
# ## Chain of layers
# ------------------------------------------------------------------------------------------

# Let's define a chain of layers
struct Chain
    layers
    Chain(layers...) = new(layers)
end
(c::Chain)(x) = (for l in c.layers; x = l(x); end; x)
(c::Chain)(x,y) = nll(c(x),y)

# ------------------------------------------------------------------------------------------
# ## Dense layers
# ------------------------------------------------------------------------------------------

# Redefine dense layer (See mlp.ipynb):
struct Dense; w; b; f; end
Dense(i::Int,o::Int,f=identity) = Dense(param(o,i), param0(o), f)
(d::Dense)(x) = d.f.(d.w * mat(x,dims=1) .+ d.b)

# ------------------------------------------------------------------------------------------
# ## Word Embeddings
# `data` has each sentence tokenized into an array of words and each word mapped to a
# `UInt16` id. To use these words as inputs to a neural network we further map each word to
# a Float32 vector. We will keep these vectors in the columns of a size (X,V) matrix where X
# is the embedding dimension and V is the vocabulary size. The vectors will be initialized
# randomly, and trained just like any other network parameter. Let's define an embedding
# layer for this purpose:
# ------------------------------------------------------------------------------------------

struct Embed; w; end
Embed(vocabsize::Int,embedsize::Int) = Embed(param(embedsize,vocabsize))
(e::Embed)(x) = e.w[:,x]

# ------------------------------------------------------------------------------------------
# This is what the words, word ids and embeddings for a sentence looks like: (note the
# identical id and embedding for the 2nd and 5th words)
# ------------------------------------------------------------------------------------------

embedlayer = Embed(length(words),8)
(w,t) = data[52855]
display(permutedims(words[w]))
display(permutedims(Int.(w)))
display(embedlayer(w))

# ------------------------------------------------------------------------------------------
# ## RNN layers
# ------------------------------------------------------------------------------------------

@doc RNN

# ------------------------------------------------------------------------------------------
# ## The three taggers: MLP, RNN, biRNN
# ------------------------------------------------------------------------------------------

# ------------------------------------------------------------------------------------------
# ## Tagger0 (MLP)
# This is what Tagger0 looks like. Every tag is predicted independently. The prediction of
# each tag only depends on the corresponding word.
#
# ------------------------------------------------------------------------------------------

# ------------------------------------------------------------------------------------------
# ## Tagger1 (RNN)
# In Tagger1, the RNN layer takes its previous output as an additional input. The prediction
# of each tag is based on words to the left.
#
# ------------------------------------------------------------------------------------------

# ------------------------------------------------------------------------------------------
# ## Tagger2 (biRNN)
# In Tagger2 there are two RNNs: the forward RNN reads the sequence from left to right, the
# backward RNN reads it from right to left. The prediction of each tag is dependent on all
# the words in the sentence.
#
# ------------------------------------------------------------------------------------------

Tagger0(vocab,embed,hidden,output)=  # MLP Tagger
    Chain(Embed(vocab,embed),Dense(embed,hidden,relu),Dense(hidden,output))
Tagger1(vocab,embed,hidden,output)=  # RNN Tagger
    Chain(Embed(vocab,embed),RNN(embed,hidden,rnnType=:relu),Dense(hidden,output))
Tagger2(vocab,embed,hidden,output)=  # biRNN Tagger
    Chain(Embed(vocab,embed),RNN(embed,hidden,rnnType=:relu,bidirectional=true),Dense(2hidden,output));

# ------------------------------------------------------------------------------------------
# ## Sequence Minibatching
# Minibatching is a bit more complicated with sequences compared to simple classification
# problems, this section can be skipped on a first reading. In addition to the input and
# minibatch sizes, there is also the time dimension to consider. To keep things simple we
# will concatenate all sentences into one big sequence, then split this sequence into equal
# sized chunks. The input to the tagger will be size (B,T) where B is the minibatch size,
# and T is the chunk size. The input to the RNN layer will be size (X,B,T) where X is the
# embedding size.
# ------------------------------------------------------------------------------------------

BATCHSIZE = 64
SEQLENGTH = 32;

function seqbatch(x,y,B,T)
    N = length(x) ÷ B
    x = permutedims(reshape(x[1:N*B],N,B))
    y = permutedims(reshape(y[1:N*B],N,B))
    d = []; for i in 0:T:N-T
        push!(d, (x[:,i+1:i+T], y[:,i+1:i+T]))
    end
    return d
end
allw = vcat((x->x[1]).(data)...)
allt = vcat((x->x[2]).(data)...)
d = seqbatch(allw, allt, BATCHSIZE, SEQLENGTH);

# ------------------------------------------------------------------------------------------
# This may be a bit more clear if we look at an example minibatch:
# ------------------------------------------------------------------------------------------

(x,y) = first(d)
words[x]

# ------------------------------------------------------------------------------------------
# ## Embedding a minibatch
# Julia indexing allows us to get the embeddings for this minibatch in one go as an (X,B,T)
# array where X is the embedding size, B is the minibatch size, and T is the subsequence
# length.
# ------------------------------------------------------------------------------------------

embedlayer = Embed(length(words),128)
summary(embedlayer(x))

# ------------------------------------------------------------------------------------------
# ## Experiments
# ------------------------------------------------------------------------------------------

# shuffle and split minibatches into train and test portions
shuffle!(d)
dtst = d[1:10]
dtrn = d[11:end]
length.((dtrn,dtst))

# For running experiments we will use the Adam algorithm which typically converges faster than SGD.
function trainresults(file,maker,savemodel)
    if (print("Train from scratch? "); readline()[1]=='y')
        model = maker()
        takeevery(n,itr) = (x for (i,x) in enumerate(itr) if i % n == 1)
        results = ((nll(model,dtst), zeroone(model,dtst))
                   for x in takeevery(100, progress(adam(model,repeat(dtrn,5)))))
        results = reshape(collect(Float32,flatten(results)),(2,:))
        Knet.save(file,"model",(savemodel ? model : nothing),"results",results)
        Knet.gc() # To save gpu memory
    else
        isfile(file) || download("http://people.csail.mit.edu/deniz/models/tutorial/$file",file)
        model,results = Knet.load(file,"model","results")
    end
    println(minimum(results,dims=2))
    return model,results
end

VOCABSIZE = length(words)
EMBEDSIZE = 128
HIDDENSIZE = 128
OUTPUTSIZE = length(tags);

# 2.35e-01  100.00%┣┫ 2780/2780 [00:13/00:13, 216.36i/s] [0.295007; 0.0972656]
t0maker() = Tagger0(VOCABSIZE,EMBEDSIZE,HIDDENSIZE,OUTPUTSIZE)
(t0,r0) = trainresults("tagger113a.jld2",t0maker,false);

# 1.49e-01  100.00%┣┫ 2780/2780 [00:19/00:19, 142.58i/s] [0.21358; 0.0616211]
t1maker() = Tagger1(VOCABSIZE,EMBEDSIZE,HIDDENSIZE,OUTPUTSIZE)
(t1,r1) = trainresults("tagger113b.jld2",t1maker,false);

# 9.37e-02  100.00%┣┫ 2780/2780 [00:25/00:25, 109.77i/s] [0.156669; 0.044043]
t2maker() = Tagger2(VOCABSIZE,EMBEDSIZE,HIDDENSIZE,OUTPUTSIZE)
(t2,r2) = trainresults("tagger113c.jld2",t2maker,true);

using Plots; default(fmt=:png,ls=:auto,ymirror=true)

plot([r0[2,:], r1[2,:], r2[2,:]]; xlabel="x100 updates", ylabel="error",
    ylim=(0,0.15), yticks=0:0.01:0.15, labels=["MLP","RNN","biRNN"])

plot([r0[1,:], r1[1,:], r2[1,:]]; xlabel="x100 updates", ylabel="loss",
    ylim=(0,.5), yticks=0:0.1:.5, labels=["MLP","RNN","biRNN"])

# ------------------------------------------------------------------------------------------
# ## Playground
# Below, you can type and tag your own sentences:
# ------------------------------------------------------------------------------------------

wdict=Dict{String,UInt16}(); for (i,w) in enumerate(words); wdict[w]=i; end
unk = UInt16(length(words))
wid(w) = get(wdict,w,unk)
function tag(tagger,s::String)
    w = permutedims(split(s))
    t = tags[(x->x[1]).(argmax(Array(tagger(wid.(w))),dims=1))]
    vcat(w,t)
end

tag(t2,readline())

max_stars_repo_name: grenkoca/JuliaTutorials

max_issues_repo_name: grenkoca/JuliaTutorials

max_forks_repo_name: grenkoca/JuliaTutorials


text:
module BorgC

const borglib = joinpath(normpath(joinpath(dirname(Base.source_path()), "..", "deps", "usr", "lib")),"borg")

typealias BORG_Problem Ptr{Void}
typealias BORG_Archive Ptr{Void}

function BORG_Problem_create(numberOfVariables::Int, numberOfObjectives::Int, numberOfConstraints::Int, fn::Function, userParams)
    c_fn = cfunction(fn, Void, (Ptr{Float64}, Ptr{Float64}, Ptr{Float64}, Ptr{Void}))

    ret = ccall(
        (:BORG_Problem_create, borglib),
        BORG_Problem,
        (Int32, Int32, Int32, Ptr{Void}, Any),
        numberOfVariables, numberOfObjectives, numberOfConstraints, c_fn, userParams)
    return ret
end

function BORG_Problem_set_bounds(problem::BORG_Problem, index::Int, lowerBound::Float64, upperBound::Float64)
    ccall(
        (:BORG_Problem_set_bounds, borglib),
        Void,
        (BORG_Problem, Int32, Float64, Float64),
        problem, index-1, lowerBound, upperBound)
end

function BORG_Problem_set_epsilon(problem::BORG_Problem, index::Int, epsilon::Float64)
    ccall(
        (:BORG_Problem_set_epsilon, "borg"),
        Void,
        (BORG_Problem, Int32, Float64),
        problem, index-1, epsilon)
end

function BORG_Algorithm_run(problem::BORG_Problem,maxEvaluations::Int)
    ret = ccall(
        (:BORG_Algorithm_run, "borg"),
        BORG_Archive,
        (BORG_Problem, Int32),
        problem, maxEvaluations)
    return ret
end

function BORG_Archive_print(archive::BORG_Archive, fp)
    ccall(
        (:BORG_Archive_print, "borg"),
        Void,
        (BORG_Archive, Ptr{Void}),
        archive, C_NULL)
end

end # module


text:
## objective
get_objective(pb::Problem) = pb.objective

"""
get_objective(pb::Problem, pt::Point)

Return the polynomial objective evaluated at `pt`.
"""
get_objective(pb::Problem, pt::Point) = evaluate(pb.objective, pt)

set_objective!(pb::Problem, p::Polynomial) = pb.objective = p

## variables
get_variables(pb::Problem) = pb.variables

function get_variabletype(pb::Problem, varName::String)
  if !haskey(pb.variables, varName)
    error("get_variable(): Current Problem has no variable named ", varName)
  end
  pb.variables[varName]
end

has_variable(pb::Problem, var::Variable) = haskey(pb.variables, var.name) && pb.variables[var.name] == var.kind

function add_variable!(pb::Problem, x::Variable)
  if haskey(pb.variables, x.name) && !(x.kind <: pb.variables[x.name])
    error("add_variable!(): Attempting to add variable ", x.name, " (", x.kind, ") when ", x, " (", get_variabletype(pb, x.name), ") already exists.")
  end
  pb.variables[x.name] = x.kind
end

function add_variable!(pb::Problem, x::Pair{String, Type})
  return add_variable!(pb, Variable(x[1], x[2]))
end

## constraints
has_constraint(pb::Problem, cstrName::String) = haskey(pb.constraints, cstrName)

get_constraints(pb::Problem) = pb.constraints

"""
  cstr = get_constraint(pb::Problem, cstrname::String)

Return the `cstrname` constraint from `pb`.
"""
get_constraint(pb::Problem, cstrname::String) = pb.constraints[cstrname]

"""
  cstr = get_constraint(pb::Problem, pt::Point)

Return the point (dict) of constraints names to their body's value when
evaluated at `pt`.
"""
function get_constraints(pb::Problem, pt::Point)
  img = Point()
  for (cstrName, cstr) in pb.constraints
    img[Variable(cstrName, Complex)] = evaluate(cstr.p, pt)
  end
  return img
end

"""
  add_constraint!(pb::Problem, cstrName::String, cstr::Constraint)

Add the constraint `cstr` under the `cstrname` name in the `pb` problem.
"""
function add_constraint!(pb::Problem, cstrName::String, cstr::Constraint)
  if haskey(pb.constraints, cstrName)
    warn("add_constraint!(): A constraint with that name already exists ($cstrName)")
  end
  if real(cstr.ub) < real(cstr.lb) || imag(cstr.ub) < imag(cstr.lb)
    warn("add_constraint!(): ", cstrName, " Lower bound is higher than upper bound ($(cstr.lb) - $(cstr.ub))")
  end
  pb.constraints[cstrName] = cstr
end

rm_constraint!(pb::Problem, cstrName::String) = pop!(pb.constraints, cstrName)

"""
  pt = get_slacks(pb::Problem, pt::Point)

Return a point associating a constraint name to its slack, i.e. the minimum
algebraic distance between the real and imaginary parts of the body's value at
`pt` and its bounds.
"""
function get_slacks(pb::Problem, pt::Point)
  var_arr = Variable[]
  val_arr = Complex[]
  for (cstrName, cstr) in pb.constraints
    val = evaluate(cstr.p, pt)
    push!(var_arr, Variable(cstrName, Complex))
    push!(val_arr, min(real(val-cstr.lb), real(cstr.ub-val)) + min(imag(val-cstr.lb), imag(cstr.ub-val))*im)
  end
  return Point(var_arr, val_arr)
end

function get_minslack(pb::Problem, pt::Point)
  minSlack = +Inf
  minCstrName = ""
  slacks = get_slacks(pb, pt)
  for (cstrName, slack) in slacks
    if minSlack > min(real(slack), imag(slack))
      minSlack = min(real(slack), imag(slack))
      minCstrName = cstrName
    end
  end
  return minSlack, minCstrName
end

max_stars_repo_name: Gimaju/GOC

max_issues_repo_name: Gimaju/GOC

max_forks_repo_name: Gimaju/GOC


text:
sandboxed_test("attr", "test parse diff options") do test_repo, path
    opts = parse_git_diff_options(AnyDict(:context_lines=>1, :interhunk_lines=>1))
    @test opts.context_lines == 1
    @test opts.interhunk_lines == 1
end

sandboxed_test("attr", "test with oid string") do test_repo, path
  d = diff(test_repo, "605812a", "370fe9ec22",
           AnyDict(:context_lines=>1, :interhunk_lines=>1))
  @test isa(d, GitDiff)

  ds = deltas(d)
  ps = patches(d)

  hs = DiffHunk[]
  for p in ps
      hks = hunks(p)
      if hks == nothing
          continue
      end
      for h in hks
          push!(hs, h)
      end
  end
  ls = vcat([lines(h) for h in hs]...)

  @test length(d) == 5
  @test length(ds) == 5
  @test length(ps) == 5

  @test sum(x -> x.status == :added? 1 : 0, ds) == 2
  @test sum(x -> x.status == :deleted? 1 : 0, ds) == 1
  @test sum(x -> x.status == :modified? 1 : 0, ds) == 2

  @test length(hs) == 5
  @test (7 + 24 + 1 + 6 + 6) == length(ls)
  @test  sum(x -> x.line_origin == :context? 1 : 0, ls) == 1
  @test  sum(x -> x.line_origin == :addition? 1 : 0, ls) == (24 + 1 + 5 + 5)
  @test  sum(x -> x.line_origin == :deletion? 1 : 0, ls) == (7 + 1)
end

#= TODO:
sandboxed_test("attr", "test delta status char") do test_repo, path
    d  = diff(test_repo, "605812a", "370fe9ec22", AnyDict(:context_lines=>1, :interhunk_lines=>1))
    ds = deltas(d)
    @test "D" == ds[1].status_char
    @test "A" == ds[2].status_char
    @test "A" == ds[3].status_char
    @test "M" == ds[4].status_char
    @test "M" == ds[5].status_char
end
=#

sandboxed_test("attr", "test with nothing on right side") do test_repo, path
  d = diff(test_repo, "605812a", nothing, AnyDict(:context_lines=>1, :interhunk_lines=>1))
  @test isa(d, GitDiff)

  ds = deltas(d)
  ps = patches(d)

  hs = DiffHunk[]
  for p in ps
      hks = hunks(p)
      if hks == nothing
          continue
      end
      for h in hks
          push!(hs, h)
      end
  end
  ls = vcat([lines(h) for h in hs]...)

  @test length(d) == 16
  @test length(ds) == 16
  @test length(ps) == 16

  @test sum(x -> x.status == :added? 1 : 0, ds) == 0
  @test sum(x -> x.status == :deleted? 1 : 0, ds) == 16
  @test sum(x -> x.status == :modified? 1 : 0, ds) == 0

  @test length(hs) == 15

  @test length(ls) == 115
  @test sum(x -> x.line_origin == :context? 1 : 0, ls) == 0
  @test sum(x -> x.line_origin == :addition? 1 : 0, ls) == 0
  @test sum(x -> x.line_origin == :deletion? 1 : 0, ls) == 113
end

sandboxed_test("attr", "test with nothing on left side") do test_repo, path
  d = diff(test_repo, nothing, "605812a", AnyDict(:context_lines=>1, :interhunk_lines=>1))
  @test isa(d, GitDiff)

  ds = deltas(d)
  ps = patches(d)

  hs = DiffHunk[]
  for p in ps
      hks = hunks(p)
      if hks == nothing
          continue
      end
      for h in hks
          push!(hs, h)
      end
  end
  ls = vcat([lines(h) for h in hs]...)

  @test length(d) == 16
  @test length(ds) == 16
  @test length(ps) == 16

  @test sum(x -> x.status == :added? 1 : 0, ds) == 16
  @test sum(x -> x.status == :deleted? 1 : 0, ds) == 0
  @test sum(x -> x.status == :modified? 1 : 0, ds) == 0

  @test length(hs) == 15

  @test length(ls) == 115
  @test sum(x -> x.line_origin == :context? 1 : 0, ls) == 0
  @test sum(x -> x.line_origin == :addition? 1 : 0, ls) == 113
  @test sum(x -> x.line_origin == :deletion? 1 : 0, ls) == 0
end

sandboxed_test("status", "test basic diff") do test_repo, path
    d = diff_workdir(
             test_repo, "26a125ee1bf",
             AnyDict(:context_lines => 3,
                  :interhunk_lines => 1,
                  :include_ignored => true,
                  :include_untracked => true))
    ds = deltas(d)
    ps = patches(d)

    hs = DiffHunk[]
    for p in ps
        hks = hunks(p)
        if hks == nothing
            continue
        end
        for h in hks
            push!(hs, h)
        end
    end
    ls = vcat([lines(h) for h in hs]...)

    @test length(ds) == 14
    @test length(ps) == 14

    @test sum(x -> x.status == :added? 1 : 0, ds) == 0
    @test sum(x -> x.status == :deleted? 1 : 0, ds) == 4
    @test sum(x -> x.status == :modified? 1 : 0, ds) == 4
    @test sum(x -> x.status == :ignored? 1 : 0, ds) == 1
    @test sum(x -> x.status == :untracked? 1 : 0, ds) == 5

    @test length(hs) == 8

    @test length(ls) == 13
    @test sum(x -> x.line_origin == :context? 1 : 0, ls) == 4
    @test sum(x -> x.line_origin == :addition? 1 : 0, ls) == 5
    @test sum(x -> x.line_origin == :deletion? 1 : 0, ls) == 4
end

sandboxed_test("attr", "test diff with parent") do test_repo, path
    c = lookup_commit(test_repo, "605812a")
    @test isa(c, GitCommit)

    d = diff(test_repo, c, AnyDict(:context_lines=>1, :interhunk_lines=>1))
    @test isa(d, GitDiff)

    ds = deltas(d)
    ps = patches(d)

    hs = DiffHunk[]
    for p in ps
        hks = hunks(p)
        if hks == nothing
            continue
        end
        for h in hks
            push!(hs, h)
        end
    end
    ls = vcat([lines(h) for h in hs]...)

    @test length(d)  == 5
    @test length(ds) == 5
    @test length(ps) == 5

    @test sum(x -> x.status == :added? 1 : 0, ds) == 0
    @test sum(x -> x.status == :deleted? 1 : 0, ds) == 3
    @test sum(x -> x.status == :modified? 1 : 0, ds) == 2

    @test length(hs) == 4

    @test length(ls) == 51
    @test sum(x -> x.line_origin == :context? 1 : 0, ls) == 2
    @test sum(x -> x.line_origin == :addition? 1 : 0, ls) == 3
    @test sum(x -> x.line_origin == :deletion? 1 : 0, ls) == 46
end

sandboxed_test("attr", "test diff with parent for init commit") do test_repo, path
    c = lookup_commit(test_repo,
                Oid("6bab5c79cd5140d0f800917f550eb2a3dc32b0da"))
    @test isa(c, GitCommit)

    d = diff(test_repo, c,
        AnyDict(:context_lines => 1, :interhunk_lines => 1, :reverse => true))
    @test isa(d, GitDiff)

    ds = deltas(d)
    ps = patches(d)

    hs = DiffHunk[]
    for p in ps
        hks = hunks(p)
        if hks == nothing
            continue
        end
        for h in hks
            push!(hs, h)
        end
    end
    ls = vcat([lines(h) for h in hs]...)

    @test length(d)  == 13
    @test length(ds) == 13
    @test length(ps) == 13

    @test sum(x -> x.status == :added? 1 : 0, ds) == 13
    @test sum(x -> x.status == :deleted? 1 : 0, ds) == 0
    @test sum(x -> x.status == :modified? 1 : 0, ds) == 0

    @test length(hs) == 13

    @test length(ls) == 72
    @test sum(x -> x.line_origin == :context? 1 : 0, ls) == 0
    @test sum(x -> x.line_origin == :addition? 1 : 0, ls) == 70
    @test sum(x -> x.line_origin == :deletion? 1 : 0, ls) == 0
end

sandboxed_test("status", "test commit to workdir basic diff") do test_repo, path
    c = lookup_commit(test_repo, "26a125ee1bf")
    @test isa(c, GitCommit)

    d= diff_workdir(test_repo, c,
     AnyDict(:context_lines => 3,
          :interhunk_lines => 1,
          :include_ignored => true,
          :include_untracked => true))

    @test isa(d, GitDiff)

    ds = deltas(d)
    ps = patches(d)

    hs = DiffHunk[]
    for p in ps
        hks = hunks(p)
        if hks == nothing
            continue
        end
        for h in hks
            push!(hs, h)
        end
    end
    ls = vcat([lines(h) for h in hs]...)

    @test length(ds) == 14
    @test length(ps) == 14

    @test sum(x -> x.status == :added? 1 : 0, ds) == 0
    @test sum(x -> x.status == :deleted? 1 : 0, ds) == 4
    @test sum(x -> x.status == :modified? 1 : 0, ds) == 4
    @test sum(x -> x.status == :ignored? 1 : 0, ds) == 1
    @test sum(x -> x.status == :untracked? 1 : 0, ds) == 5

    @test length(hs) == 8

    @test length(ls) == 13
    @test sum(x -> x.line_origin == :context? 1 : 0, ls) == 4
    @test sum(x -> x.line_origin == :addition? 1 : 0, ls) == 5
    @test sum(x -> x.line_origin == :deletion? 1 : 0, ls) == 4
end

sandboxed_test("status", "test tree to workdir basic diff") do test_repo, path
    c = lookup_commit(test_repo, "26a125ee1bf")
    t = GitTree(c)
    @test isa(t, GitTree)

    d = diff_workdir(test_repo, t,
     AnyDict(:context_lines => 3,
             :interhunk_lines => 1,
             :include_ignored => true,
             :include_untracked => true))
    @test isa(d, GitDiff)

    ds = deltas(d)
    ps = patches(d)

    hs = DiffHunk[]
    for p in ps
        hks = hunks(p)
        if hks == nothing
            continue
        end
        for h in hks
            push!(hs, h)
        end
    end
    ls = vcat([lines(h) for h in hs]...)

    @test length(ds) == 14
    @test length(ps) == 14

    @test sum(x -> x.status == :added? 1 : 0, ds) == 0
    @test sum(x -> x.status == :deleted? 1 : 0, ds) == 4
    @test sum(x -> x.status == :modified? 1 : 0, ds) == 4
    @test sum(x -> x.status == :ignored? 1 : 0, ds) == 1
    @test sum(x -> x.status == :untracked? 1 : 0, ds) == 5

    @test length(hs) == 8

    @test length(ls) == 13
    @test sum(x -> x.line_origin == :context? 1 : 0, ls) == 4
    @test sum(x -> x.line_origin == :addition? 1 : 0, ls) == 5
    @test sum(x -> x.line_origin == :deletion? 1 : 0, ls) == 4
end

sandboxed_test("status", "test tree to workdir diff merge") do test_repo, path
    idx = GitIndex(test_repo)
    c = lookup_commit(test_repo, "26a125ee1bf")
    t = GitTree(c)

    # merge diffs to simulate "git diff 26a125ee1bf"
    diff1 = diff(test_repo, t, idx,
                AnyDict(:include_ignored=>true, :include_untracked=>true))
    diff2 = diff(test_repo, idx,
                AnyDict(:include_ignored=>true, :include_untracked=>true))
    merge!(diff1, diff2)

    ds = deltas(diff1)
    ps = patches(diff1)

    hs = DiffHunk[]
    for p in ps
        hks = hunks(p)
        if hks == nothing
            continue
        end
        for h in hks
            push!(hs, h)
        end
    end
    ls = vcat([lines(h) for h in hs]...)

    @test length(ds) == 15
    @test length(ps) == 15

    @test sum(x -> x.status == :added? 1 : 0, ds) == 2
    @test sum(x -> x.status == :deleted? 1 : 0, ds) == 5
    @test sum(x -> x.status == :modified? 1 : 0, ds) == 4
    @test sum(x -> x.status == :ignored? 1 : 0, ds) == 1
    @test sum(x -> x.status == :untracked? 1 : 0, ds) == 3

    @test length(hs) == 11

    @test length(ls) == 17
    @test sum(x -> x.line_origin == :context? 1 : 0, ls) == 4
    @test sum(x -> x.line_origin == :addition? 1 : 0, ls) == 8
    @test sum(x -> x.line_origin == :deletion? 1 : 0, ls) == 5
end

sandboxed_test("status", "test tree to workdir diff stats") do test_repo, path
    idx = GitIndex(test_repo)
    c = lookup_commit(test_repo, "26a125ee1bf")
    t = GitTree(c)

    # merge diffs to simulate "git diff 26a125ee1bf"
    diff1 = diff(test_repo, t, idx,
                AnyDict(:include_ignored=>true, :include_untracked=>true))
    diff2 = diff(test_repo, idx,
                AnyDict(:include_ignored=>true, :include_untracked=>true))
    merge!(diff1, diff2)
    diff_stat = stat(diff1)

    @test diff_stat.files == 11
    @test diff_stat.adds  == 8
    @test diff_stat.dels  == 5

    expected_patch_stat = Any[
      [ 0, 1, 1 ], [ 1, 0, 2 ], [ 1, 0, 2 ], [ 0, 1, 1 ], [ 2, 0, 3 ],
      [ 0, 1, 1 ], [ 0, 1, 1 ], [ 1, 0, 1 ], [ 2, 0, 2 ], [ 0, 1, 1 ],
      [ 1, 0, 2 ]
    ]

    i = 1
    for patch in patches(diff1)
        s = delta(patch).status
        if s == :unmodified || s == :ignored || s == :untracked
            continue
        end

        @test isa(patch, GitPatch)
        expected_adds  = expected_patch_stat[i][1]
        expected_dels  = expected_patch_stat[i][2]
        expected_lines = expected_patch_stat[i][3]

        patch_stat = stat(patch)

        @test expected_adds == patch_stat.adds
        @test expected_dels == patch_stat.dels
        @test nchanges(patch) == patch_stat.adds + patch_stat.dels
        #@test expected_lines == nlines(patch)
        i += 1
    end
end

sandboxed_test("attr", "test tree to tree basic diff") do test_repo, path
    a = GitTree(lookup_commit(test_repo, "605812a"))
    b = GitTree(lookup_commit(test_repo, "370fe9ec22"))
    c = GitTree(lookup_commit(test_repo, "f5b0af1fb4f5c"))

    d = diff(test_repo, a, b, AnyDict(:context_lines => 1, :interhunk_lines => 1))
    ds = deltas(d)
    ps = patches(d)
    hs = DiffHunk[]
    for p in ps
        hks = hunks(p)
        if hks == nothing
            continue
        end
        for h in hks
            push!(hs, h)
        end
    end
    ls = vcat([lines(h) for h in hs]...)

    @test length(d) == 5
    @test length(ds) == 5
    @test length(ps) == 5

    @test sum(x -> x.status == :added? 1 : 0, ds) == 2
    @test sum(x -> x.status == :deleted? 1 : 0, ds) == 1
    @test sum(x -> x.status == :modified? 1 : 0, ds) == 2

    @test length(hs) == 5

    @test length(ls) == (7 + 24 + 1 + 6 + 6)
    @test sum(x -> x.line_origin == :context? 1 : 0, ls) == 1
    @test sum(x -> x.line_origin == :addition? 1 : 0, ls) == (24 + 1 + 5 + 5)
    @test sum(x -> x.line_origin == :deletion? 1 : 0, ls) == (7 + 1)

    d = diff(test_repo, c, b, AnyDict(:context_lines => 1, :interhunk_lines => 1))
    ds = deltas(d)
    ps = patches(d)
    hs = DiffHunk[]
    for p in ps
        hks = hunks(p)
        if hks == nothing
            continue
        end
        for h in hks
            push!(hs, h)
        end
    end
    ls = vcat([lines(h) for h in hs]...)

    @test length(ds) == 2
    @test length(ps) == 2

    @test sum(x -> x.status == :added? 1 : 0, ds) == 0
    @test sum(x -> x.status == :deleted? 1 : 0, ds) == 0
    @test sum(x -> x.status == :modified? 1 : 0, ds) == 2

    @test length(hs) == 2

    @test length(ls) == (8 + 15)
    @test sum(x -> x.line_origin == :context? 1 : 0, ls) == 1
    @test sum(x -> x.line_origin == :addition? 1 : 0, ls) == 1
    @test sum(x -> x.line_origin == :deletion? 1 : 0, ls) == (7 + 14)
end

sandboxed_test("attr", "test tree to tree with empty tree") do test_repo, path
    a = GitTree(lookup_commit(test_repo, "605812a"))
    d = diff(test_repo, a, nothing, AnyDict(:context_lines => 1, :interhunk_lines => 1))
    ds = deltas(d)
    ps = patches(d)
    hs = DiffHunk[]
    for p in ps
        hks = hunks(p)
        if hks == nothing
            continue
        end
        for h in hks
            push!(hs, h)
        end
    end
    ls = vcat([lines(h) for h in hs]...)

    @test length(d) == 16
    @test length(ds) == 16
    @test length(ps) == 16

    @test sum(x -> x.status == :added? 1 : 0, ds) == 0
    @test sum(x -> x.status == :deleted? 1 : 0, ds) == 16
    @test sum(x -> x.status == :modified? 1 : 0, ds) == 0

    @test length(hs) == 15

    @test length(ls) == 115
    @test sum(x -> x.line_origin == :context? 1 : 0, ls) == 0
    @test sum(x -> x.line_origin == :addition? 1 : 0, ls) == 0
    @test sum(x -> x.line_origin == :deletion? 1 : 0, ls) == 113
end

sandboxed_test("attr", "test tree to tree with rev string") do test_repo, path
    a = GitTree(lookup_commit(test_repo, "605812a"))
    d = diff(test_repo, a, "370fe9ec22",
                AnyDict(:context_lines => 1, :interhunk_lines => 1))
    @test isa(d, GitDiff)
    ds = deltas(d)
    ps = patches(d)
    hs = DiffHunk[]
    for p in ps
        hks = hunks(p)
        if hks == nothing
            continue
        end
        for h in hks
            push!(hs, h)
        end
    end
    ls = vcat([lines(h) for h in hs]...)

    @test length(d) == 5
    @test length(ds) == 5
    @test length(ps) == 5

    @test sum(x -> x.status == :added? 1 : 0, ds) == 2
    @test sum(x -> x.status == :deleted? 1 : 0, ds) == 1
    @test sum(x -> x.status == :modified? 1 : 0, ds) == 2

    @test length(hs) == 5

    @test length(ls) == (7 + 24 + 1 + 6 + 6)
    @test sum(x -> x.line_origin == :context? 1 : 0, ls) == 1
    @test sum(x -> x.line_origin == :addition? 1 : 0, ls) == (24 + 1 + 5 + 5)
    @test sum(x -> x.line_origin == :deletion? 1 : 0, ls) == (7 + 1)
end

sandboxed_test("attr", "test tree to tree with diff merge") do test_repo, path
    a = GitTree(lookup_commit(test_repo, "605812a"))
    b = GitTree(lookup_commit(test_repo, "370fe9ec22"))
    c = GitTree(lookup_commit(test_repo, "f5b0af1fb4f5c"))

    ab = diff(test_repo, a, b)
    cb = diff(test_repo, c, b)
    d  = merge!(ab, cb)

    ds = deltas(d)
    ps = patches(d)
    hs = DiffHunk[]
    for p in ps
        hks = hunks(p)
        if hks == nothing
            continue
        end
        for h in hks
            push!(hs, h)
        end
    end
    ls = vcat([lines(h) for h in hs]...)

    @test length(d) == 6
    @test length(ds) == 6
    @test length(ps) == 6

    @test sum(x -> x.status == :added? 1 : 0, ds) == 2
    @test sum(x -> x.status == :deleted? 1 : 0, ds) == 1
    @test sum(x -> x.status == :modified? 1 : 0, ds) == 3

    @test length(hs) == 6

    @test length(ls) == 59
    @test sum(x -> x.line_origin == :context? 1 : 0, ls) == 1
    @test sum(x -> x.line_origin == :addition? 1 : 0, ls) == 36
    @test sum(x -> x.line_origin == :deletion? 1 : 0, ls) == 22
end

sandboxed_test("unsymlinked.git", "test symlink blob mode changed to reg file") do test_repo, path
    a = GitTree(lookup_commit(test_repo, "7fccd7"))
    b = GitTree(lookup_commit(test_repo, "806999"))

    d = diff(test_repo, a, b)

    ds = deltas(d)
    ps = patches(d)

    @test length(d) == 3
    @test length(ds) == 3
    @test length(ps) == 3

    @test sum(x -> x.status == :added? 1 : 0, ds) == 1
    @test sum(x -> x.status == :deleted? 1 : 0, ds) == 2
    @test sum(x -> x.status == :modified? 1 : 0, ds) == 0
    @test sum(x -> x.status == :typechange? 1 : 0, ds) == 0
end

sandboxed_test("unsymlinked.git", "test symlink blob mode changed to reg file as typechange") do test_repo, path
    a = GitTree(lookup_commit(test_repo, "7fccd7"))
    b = GitTree(lookup_commit(test_repo, "806999"))

    d = diff(test_repo, a, b, AnyDict(:include_typechange => true))

    ds = deltas(d)
    ps = patches(d)

    @test length(d) == 2
    @test length(ds) == 2
    @test length(ps) == 2

    @test sum(x -> x.status == :added? 1 : 0, ds) == 0
    @test sum(x -> x.status == :deleted? 1 : 0, ds) == 1
    @test sum(x -> x.status == :modified? 1 : 0, ds) == 0
    @test sum(x -> x.status == :typechange? 1 : 0, ds) == 1
end

sandboxed_test("unsymlinked.git", "test regular blob mode change to exe file") do test_repo, path
    a = GitTree(lookup_commit(test_repo, "806999"))
    b = GitTree(lookup_commit(test_repo, "a8595c"))

    d = diff(test_repo, a, b)

    ds = deltas(d)
    ps = patches(d)

    @test length(d) == 1
    @test length(ds) == 1
    @test length(ps) == 1

    @test sum(x -> x.status == :added? 1 : 0, ds) == 0
    @test sum(x -> x.status == :deleted? 1 : 0, ds) == 0
    @test sum(x -> x.status == :modified? 1 : 0, ds) == 1
    @test sum(x -> x.status == :typechange? 1 : 0, ds) == 0
end

sandboxed_test("diff", "test size") do test_repo, path
    a = lookup(test_repo, Oid("d70d245ed97ed2aa596dd1af6536e4bfdb047b69"))
    b = lookup(test_repo, Oid("7a9e0b02e63179929fed24f0a3e0f19168114d10"))

    d = diff(test_repo, GitTree(a), GitTree(b))
    @test length(d) == 2
end

sandboxed_test("diff", "test each delta") do test_repo, path
    a = lookup(test_repo, Oid("d70d245ed97ed2aa596dd1af6536e4bfdb047b69"))
    b = lookup(test_repo, Oid("7a9e0b02e63179929fed24f0a3e0f19168114d10"))

    d = diff(test_repo, GitTree(a), GitTree(b))
    ds = deltas(d)

    @test length(d) == 2

    @test ds[1].old_file.path == "another.txt"
    @test ds[1].new_file.path == "another.txt"

    @test ds[1].isbinary == false

    @test ds[2].old_file.path == "readme.txt"
    @test ds[2].new_file.path == "readme.txt"

    @test ds[2].isbinary == false
end

sandboxed_test("diff", "test diff treats files bigger than max size as binary") do test_repo, path
    a = lookup(test_repo, Oid("d70d245ed97ed2aa596dd1af6536e4bfdb047b69"))
    b = lookup(test_repo, Oid("7a9e0b02e63179929fed24f0a3e0f19168114d10"))
    d = diff(test_repo, GitTree(a), GitTree(b), AnyDict(:max_size=>10))
    @test length(patches(d)) == 2
    @test patch(d) == replace("diff --git a/another.txt b/another.txt
index 3e5bcba..546c735 100644
Binary files a/another.txt and b/another.txt differ
diff --git a/readme.txt b/readme.txt
index 7b808f7..29ab705 100644
Binary files a/readme.txt and b/readme.txt differ
", "\r\n", "\n")
end

sandboxed_test("diff", "test constraining paths") do test_repo, path
    a = GitTree(lookup(test_repo, Oid("d70d245ed97ed2aa596dd1af6536e4bfdb047b69")))
    b = GitTree(lookup(test_repo, Oid("7a9e0b02e63179929fed24f0a3e0f19168114d10")))

    d = diff(test_repo, a, b, AnyDict(:paths => ["readme.txt"]))
    @test "M\treadme.txt\n" == patch(d, format=:name_status)

    d = diff(test_repo, a, b, AnyDict(:paths => ["r*.txt"]))
    @test "M\treadme.txt\n" == patch(d, format=:name_status)

    d = diff(test_repo, a, b, AnyDict(:paths => ["*.txt"]))
    @test "M\tanother.txt\nM\treadme.txt\n" == patch(d, format=:name_status)

    d = diff(test_repo, a, b, AnyDict(:paths => ["*.txt"], :disable_pathspec_match => true))
    @test "" == patch(d, format=:name_status)

    d = diff(test_repo, a, b, AnyDict(:paths => ["readme.txt"], :disable_pathspec_match => true))
    @test "M\treadme.txt\n" == patch(d, format=:name_status)
end

sandboxed_test("diff", "test patch") do test_repo, path
    a = GitTree(lookup(test_repo, Oid("d70d245ed97ed2aa596dd1af6536e4bfdb047b69")))
    b = GitTree(lookup(test_repo, Oid("7a9e0b02e63179929fed24f0a3e0f19168114d10")))
    d = diff(test_repo, a, b, AnyDict(:context_lines => 0))

    @test patch(d) == replace("diff --git a/another.txt b/another.txt
index 3e5bcba..546c735 100644
--- a/another.txt
+++ b/another.txt
@@ -2 +2 @@ Git is fast. With Git, nearly all operations are performed locally, giving
-it a huge speed advantage on centralized systems that constantly have to
+it an huge speed advantage on centralized systems that constantly have to
@@ -11,4 +10,0 @@ from the start.
-Let's see how common operations stack up against Subversion, a common
-centralized version control system that is similar to CVS or
-Perforce. Smaller is faster.
-
@@ -34,0 +31,4 @@ SVN.
+Let's see how common operations stack up against Subversion, a common
+centralized version control system that is similar to CVS or
+Perforce. Smaller is faster.
+
diff --git a/readme.txt b/readme.txt
index 7b808f7..29ab705 100644
--- a/readme.txt
+++ b/readme.txt
@@ -1 +1 @@
-The Git feature that really makes it stand apart from nearly every other SCM
+The Git feature that r3ally mak3s it stand apart from n3arly 3v3ry other SCM
@@ -10,4 +9,0 @@ This means that you can do things like:
-Frictionless Context Switching. Create a branch to try out an idea, commit a
-few times, switch back to where you branched from, apply a patch, switch
-back to where you are experimenting, and merge it in.
-
@@ -27,3 +22,0 @@ Notably, when you push to a remote repository, you do not have to push all
-of your branches. You can choose to share just one of your branches, a few
-of them, or all of them. This tends to free people to try new ideas without
-worrying about having to plan how and when they are going to merge it in or
@@ -35 +28 @@ incredibly easy and it changes the way most developers work when they learn
-it.
+it.!
\\ No newline at end of file
", "\r\n", "\n")
end

sandboxed_test("diff", "test iteration") do test_repo, path
    a = GitTree(lookup(test_repo, Oid("d70d245ed97ed2aa596dd1af6536e4bfdb047b69")))
    b = GitTree(lookup(test_repo, Oid("7a9e0b02e63179929fed24f0a3e0f19168114d10")))
    d = diff(test_repo, a, b)

    ps = patches(d)
    @test length(ps) == 2

    @test "another.txt" == delta(ps[1]).old_file.path
    @test "another.txt" == delta(ps[1]).new_file.path

    @test delta(ps[1]).isbinary == false

    hs = hunks(ps[1])
    @test length(hs) == 3

    @test beginswith(hs[1].header, "@@ -1,5 +1,5 @@")
    @test beginswith(hs[2].header, "@@ -8,10 +8,6 @@")
    @test beginswith(hs[3].header, "@@ -32,6 +28,10 @@")

    @test "readme.txt" == delta(ps[2]).old_file.path
    @test "readme.txt" == delta(ps[2]).new_file.path

    @test delta(ps[2]).isbinary == false

    hs = hunks(ps[2])
    @test length(hs) == 3

    @test beginswith(hs[1].header, "@@ -1,4 +1,4 @@")
    @test beginswith(hs[2].header, "@@ -7,10 +7,6 @@")
    @test beginswith(hs[3].header, "@@ -24,12 +20,9 @@")

    ls = lines(hs[1])
    @test length(ls) == 5

    @test :deletion == ls[1].line_origin
    @test "The Git feature that really makes it stand apart from nearly every other SCM\n" == ls[1].content
    @test 0 == ls[1].content_offset

    @test :addition == ls[2].line_origin
    @test "The Git feature that r3ally mak3s it stand apart from n3arly 3v3ry other SCM\n" == ls[2].content
    @test 0 == ls[2].content_offset

    @test :context == ls[3].line_origin
    @test "out there is its branching model.\n" == ls[3].content
    @test ls[3].content_offset == nothing

    @test :context == ls[4].line_origin
    @test "\n" == ls[4].content

    @test :context == ls[5].line_origin
    @test "Git allows and encourages you to have multiple local branches that can be\n" == ls[5].content

    ls = lines(hs[3])

    @test length(ls) == 14

    @test :deletion == ls[4].line_origin
    @test "of your branches. You can choose to share just one of your branches, a few\n" == ls[4].content
    @test 1248 == ls[4].content_offset

    @test :deletion == ls[5].line_origin
    @test "of them, or all of them. This tends to free people to try new ideas without\n" == ls[5].content
    @test 1323 == ls[5].content_offset

    @test :deletion == ls[12].line_origin
    @test "it.\n" == ls[12].content
    @test 1721 == ls[12].content_offset

    @test :addition == ls[13].line_origin
    @test "it.!" == ls[13].content
    @test 1289 == ls[13].content_offset
end

sandboxed_test("diff", "test patch compact") do test_repo, path
    a = GitTree(test_repo[Oid("d70d245ed97ed2aa596dd1af6536e4bfdb047b69")])
    b = GitTree(test_repo[Oid("7a9e0b02e63179929fed24f0a3e0f19168114d10")])
    d = diff(test_repo, a, b)
    @test patch(d, format=:name_status) == replace("M\tanother.txt
M\treadme.txt
", "\r\n", "\n")
end

sandboxed_test("diff", "test stats") do test_repo, path
    a = GitTree(test_repo[Oid("d70d245ed97ed2aa596dd1af6536e4bfdb047b69")])
    b = GitTree(test_repo[Oid("7a9e0b02e63179929fed24f0a3e0f19168114d10")])
    d = diff(test_repo, a, b)

    diff_stat = stat(d)
    @test diff_stat.files == 2
    @test diff_stat.adds  == 7
    @test diff_stat.dels  == 14

    expected_patch_stat = Any[[ 5, 5, 26 ], [ 2, 9, 28 ]]

    for (patch, expected) in zip(patches(d), expected_patch_stat)
        patch_stat = stat(patch)
        @test expected[1] == patch_stat.adds
        @test expected[2] == patch_stat.dels
        @test expected[3] == nlines(patch)
    end
end


text:
## This implements data types and methods for direct detection imaging (i.e. standard array imaging at image plane)

"""
    ArrayDD <: BAMFData

data type for direct detection (i.e. a single image at a camera).
"""
abstract type ArrayDD <: BAMFData end

mutable struct ArrayDD_CUDA <: ArrayDD  # direct detection data
    sz::Int
    data::CuArray{Float32,2}
end
mutable struct ArrayDD_CPU <: ArrayDD  # direct detection data
    sz::Int
    data::Array{Float32,2}
end

ArrayDD_CUDA(sz) = ArrayDD(sz, CuArray{Float32}(undef, sz, sz))
ArrayDD(sz) = ArrayDD_CPU(sz, Array{Float32}(undef, sz, sz))
# function deepcopy(a::ArrayDD_CPU)
#     b=ArrayDD(a.sz)
#     for nn=1:b.sz*b.sz
#         b.data[nn]=a.data[nn]
#     end
# end

function genmodel!(m::StateFlatBg,rjs::RJStruct,model::ArrayDD)
    genmodel!(m, model.sz, rjs.psf, model.data)
end

function genmodel!(m::StateFlatBg,psf::PSF,model::ArrayDD)
    genmodel!(m, model.sz, psf, model.data)
end

function genmodel!(m::StateFlatBg,psf::MicroscopePSFs.PSF,model::ArrayDD)
    genmodel!(m, model.sz, psf, model.data)
end

function likelihoodratio(m::ArrayDD, mtest::ArrayDD, d::ArrayDD)
    return likelihoodratio(m.data,mtest.data,d.data)
end

function poissrnd(d::ArrayDD)
    out = ArrayDD(d.sz)
    for nn = 1:d.sz^2
        out.data[nn] = Float32(rand(Poisson(Float64(d.data[nn]))))
    end
    return out
end

function genBAMFData(data::ArrayDD)
    return ArrayDD(data.sz)
end

function calcresiduum(model::ArrayDD,data::ArrayDD)
    residuum = ArrayDD(model.sz);
    for nn=1:model.sz*model.sz
        residuum.data[nn]=data.data[nn]-model.data[nn]
    end
    return residuum
end

"""
    makepdf!(a::ArrayDD)

make all elements sum to 1
"""
function makepdf!(a::ArrayDD)
mysum=0;
for nn=1:a.sz*a.sz
    a.data[nn]=max(0,a.data[nn])
    mysum+=a.data[nn]
end
for nn=1:a.sz*a.sz
    a.data[nn]/=mysum
end
end

"""
    makecdf!(a::ArrayData)

make the array a normalized CDF
"""
function makecdf!(a::ArrayDD)
makepdf!(a)
for nn=2:a.sz*a.sz
    a.data[nn]+=a.data[nn-1];
end
end

"""
    arrayrand!(a::ArrayDD)

Pulls random number from pdf array. This converts input to cdf
"""
function arrayrand!(a::ArrayDD)
#pull random number from pdf array
#this converts input to cdf
    makecdf!(a)
    r=rand()
    #r must be less than the max value in a.data or else the while loop reads over the array
    r=min(r,maximum(a.data))
    nn=1
    while a.data[nn]<r
        nn+=1
    end
    ii=Int(rem(nn,a.sz))
    jj=Int(ceil(nn/a.sz))
    return ii,jj
end

function arraypdf(a::ArrayDD,ii::Int,jj::Int)
#calculate probability at pixel index
mysum=0;
for nn=1:a.sz*a.sz
    mysum+=max(1,a.data[nn]);
end
  return max(1,a.data[ii,jj])/mysum
end


text:
using Test
using MPI
using ClimateMachine.Mesh.Topologies

function main()
    MPI.Init()
    comm = MPI.COMM_WORLD
    crank = MPI.Comm_rank(comm)
    csize = MPI.Comm_size(comm)

    @assert csize == 2

    FT = Float64
    Nx = 3
    Ny = 2
    x = range(FT(0); length = Nx + 1, stop = 1)
    y = range(FT(0); length = Ny + 1, stop = 1)

    topology = BrickTopology(
        comm,
        (x, y);
        boundary = ((1, 2), (3, 4)),
        periodicity = (true, true),
        connectivity = :face,
    )

    elems = topology.elems
    realelems = topology.realelems
    ghostelems = topology.ghostelems
    sendelems = topology.sendelems
    elemtocoord = topology.elemtocoord
    elemtoelem = topology.elemtoelem
    elemtoface = topology.elemtoface
    elemtoordr = topology.elemtoordr
    elemtobndy = topology.elemtobndy
    nabrtorank = topology.nabrtorank
    nabrtorecv = topology.nabrtorecv
    nabrtosend = topology.nabrtosend

    globalelemtoface = [
        2 2 2 2 2 2
        1 1 1 1 1 1
        4 4 4 4 4 4
        3 3 3 3 3 3
    ]

    globalelemtoordr = ones(Int, size(globalelemtoface))
    globalelemtobndy = zeros(Int, size(globalelemtoface))

    if crank == 0
        nrealelem = 3
        globalelems = [1, 2, 3, 4, 5, 6]
        elemtoelem_expect = [
            6 4 2 4 5 6
            5 3 4 4 5 6
            2 1 5 4 5 6
            2 1 5 4 5 6
        ]
        nabrtorank_expect = [1]
        nabrtorecv_expect = UnitRange{Int}[1:3]
        nabrtosend_expect = UnitRange{Int}[1:3]
    elseif crank == 1
        nrealelem = 3
        globalelems = [4, 5, 6, 1, 2, 3]
        elemtoelem_expect = [
            6 4 2 4 5 6
            5 3 4 4 5 6
            3 6 1 4 5 6
            3 6 1 4 5 6
        ]
        nabrtorank_expect = [0]
        nabrtorecv_expect = UnitRange{Int}[1:3]
        nabrtosend_expect = UnitRange{Int}[1:3]
    end

    @test elems == 1:length(globalelems)
    @test realelems == 1:nrealelem
    @test ghostelems == (nrealelem + 1):length(globalelems)
    @test elemtoface[:, realelems] ==
          globalelemtoface[:, globalelems[realelems]]
    @test elemtoelem == elemtoelem_expect
    @test elemtobndy == globalelemtobndy[:, globalelems]
    @test elemtoordr == ones(eltype(elemtoordr), size(elemtoordr))
    @test nabrtorank == nabrtorank_expect
    @test nabrtorecv == nabrtorecv_expect
    @test nabrtosend == nabrtosend_expect

    @test collect(realelems) ==
          sort(union(topology.exteriorelems, topology.interiorelems))
    @test unique(sort(sendelems)) == topology.exteriorelems
    @test length(intersect(topology.exteriorelems, topology.interiorelems)) == 0
end

main()


text:
function init_window(engine::Engine, def::Dict{Symbol, Any})
	GLFW.Init()
	GLFW.WindowHint(GLFW.DEPTH_BITS, get(def, :depth_bits, 24))
	GLFW.WindowHint(GLFW.STENCIL_BITS, get(def, :stencil_bits, 8))
	local ctx, api
	if engine.api == :gl
		api = GLFW.OPENGL_API
		ctx = GLFW.NATIVE_CONTEXT_API
	elseif engine.api == :gles
		api = GLFW.OPENGL_ES_API
		ctx = GLFW.EGL_CONTEXT_API
	else
		api = GLFW.NO_API
	end
	GLFW.WindowHint(GLFW.CLIENT_API, api)
	GLFW.WindowHint(GLFW.CONTEXT_CREATION_API, ctx)
	ver = map(Int, get(def, :api_version, engine.apiVersion))
	GLFW.WindowHint(GLFW.CONTEXT_VERSION_MAJOR, ver[1])
	GLFW.WindowHint(GLFW.CONTEXT_VERSION_MINOR, ver[2])
	winSize = map(Int, get(def, :window_size, [640, 480]))
	engine.window = GLFW.CreateWindow(winSize[1], winSize[2], get(def, :window_title, "GamEn.jl"))
	GLFW.MakeContextCurrent(engine.window)
	GLFW.SwapInterval(get(def, :swap_interval, 0))
end

function done_window(engine::Engine)
	GLFW.Terminate()
end

function init_viewport(engine::Engine)
	set_viewport(engine, GLFW.GetFramebufferSize(engine.window)...)
	GLFW.SetFramebufferSizeCallback(engine.window, (win::GLFW.Window, width::Cint, height::Cint) -> set_viewport(engine, width, height))
end

function set_viewport(engine::Engine, width::Integer, height::Integer)
	ModernGL.glViewport(0, 0, width, height)

	call_event(engine, :viewport, width, height)
end


text:
## Quantile Functions
## ==================

"Computes quantile of interval, quantile is monotonic function so its easy"
function quantile{T}(X::Distributions.Distribution, p::Interval{T})
  # Quantiles must be between 0 and 1
  p_valid = intersect(p, unit(p))
  lb = quantile(X, p_valid.l)
  ub = quantile(X, p_valid.u)
  Interval{T}(lb, ub)
end

# quantile(X::Distributions.Distribution, i::Id, ω::Omega) = quantile(X, ω[i])

"Quantile when μ, σ and p are all intervals. This may not be correct."
function quantile{T}(::Type{Normal}, μ::Interval{T}, σ::Interval{T}, p::Interval{T})
  # FIXME, Float specific and feels a bit hacky.
  σ_valid = intersect(σ, Interval{T}(nextfloat(zero(T)), inf(T)))
  p_valid = intersect(p, unit(p))
  lb = quantile(Normal(μ.l, σ_valid.u), p_valid.l)
  ub = quantile(Normal(μ.u, σ_valid.l), p_valid.u)
  Interval{T}(lb, ub)
end

# Lebesgue Measure
measure(i::Interval) = i.u - i.l


text:
#julia
using Makie, GeometryTypes
scene = Scene()

large_sphere = HyperSphere(Point3f0(0), 1f0)
positions = decompose(Point3f0, large_sphere)
linepos = view(positions, rand(1:length(positions), 1000))

lines(linepos, linewidth = 0.1, color = :black)
scatter(positions, strokewidth = 0.01, strokecolor = :white, color = RGBA(0.9, 0.2, 0.4, 0.4))
r = linspace(-0.1, 1.1, 5)
axis(r, r, r)
center!(scene)


text:
using TensorDSL
using Test

V3 = VectorSpace(3)
V2 = VectorSpace(2)
Vi = VectorSpace()
Vj = VectorSpace()

x = FreeIndex(V3, "x")
y = FreeIndex(V2, "y")
z = FreeIndex(Vi, "z")
w = FreeIndex(Vj, "w")
fixedindices = Indices(FixedIndex(V3, 2), FixedIndex(V2, 1))

A = VariableTensor("A", V3, V2)
B = Tensor(fill(1.2, (3, 2)), V3, V2)
C = VariableTensor("C", Vj, Vi, Vi')
D = VariableTensor("D", V2', Vi)
E = VariableTensor("E", V2', V3', Vi)
F = Tensor(fill(1.5, (2,3)), V2', V3)
Z = ZeroTensor(V3, V2)
a = VariableTensor("a")
b = VariableTensor("b")
arrayG = [cos(a), 4*sin(a), 4, 7im]
G = Tensor(arrayG, RnSpace(4))
H = VariableTensor("H", V2', RnSpace(5))
# TODO make the addition nicer
I = Tensor(reshape([1, 2, 3, 4, 5, 6, 7, 8, 9, 10, 11, ConstantTensor(11) + a], (2,2,3)), V2, V2', V3)

@testset "TensorDSL.jl" begin
    @testset "Operations" begin
        @testset "Addition" begin
            @test (A + B).shape == A.shape
            @test A + Z == A
            @test (A + Z)[1, 2] == A[1, 2]
            @test Z + Z + Z == Z
            @test (A + B).children == (A, B)
            @test (A[x,1] + B[1, y]).freeindices == (x, y)
            @test_throws DimensionMismatch A + D
            @test (A[x, y] + D[y', z]).freeindices == (x,z)
        end
        @testset "Index" begin
            @test A[x, 1].shape == ()
            @test A[x, 1].freeindices == (x,)
            @test_throws DimensionMismatch A[x, x]
            @test A[x].shape == (V2,)
            @test A[x] isa ComponentTensorOperation
            @test E[y'].shape == (V3', Vi)
        end
        @testset "Outer product" begin
            @test (A⊗B).shape == (V3, V2, V3, V2)
            @test (A⊗B).children == (A, B)
            @test (A[x,1] ⊗ B[1, y]).freeindices == (x, y)
            @test (A - B).shape == (V3, V2)
            @test (-B).shape == (V3, V2)
            @test (A - B) isa AddOperation
        end
        @testset "Power" begin
            @test A[1,2]^2 isa PowerOperation
            @test sqrt(a) isa PowerOperation
            @test sqrt(a).children == (a, ConstantTensor(1//2))
        end
        @testset "Division" begin
            @test_throws DivideError A / ZeroTensor()
            @test_throws DivideError A / 0
            @test_throws Exception A / A
            @test A / 5 isa DivisionOperation
            @test (A / 4).children == (A, ConstantTensor(4))
            @test A / 1 == A
        end
        @testset "Component tensor" begin
            @test componenttensor(A[x, 1], x).shape == (V3,)
            @test componenttensor(A[x, y], y).shape == (V2,)
            @test componenttensor(A[x, y], y, x).shape == (V2, V3)
            @test componenttensor(E[y', x', z], x', z).shape == (V3', Vi)
            @test componenttensor(E[y', x', z], x', z).freeindices == (y',)
            @test_throws DomainError componenttensor(A[x, y], x, z)
        end
        @testset "Tensor contraction" begin
            @test (A[x, y]*D[y', z]).shape == ()
            @test (A[x, y]*D[y', z]).freeindices == (x, z)
            @test (A[x, y]*E[y', x', z]).shape == ()
            @test (A[x, y]*E[y', x', z]).freeindices == (z,)
            @test componenttensor(C[w, z, z'], w).shape == (Vj,)
            @test C[w, z, z'].freeindices == (w,)
            @test C[w, z, z'].children[2] == Indices(z)
        end
        @testset "Transpose" begin
            @test transpose(A).shape == reverse(A.shape)
            @test transpose(A) isa ComponentTensorOperation
            @test transpose(C).shape == reverse(C.shape)
        end
        @testset "Trace" begin
            @test trace(VariableTensor("X", V3, V3')) isa IndexSumOperation
            @test_throws DomainError trace(VariableTensor("X", V3, V3))
            @test_throws MethodError trace(VariableTensor("X", V3, V3, V3'))
        end
    end
    @testset "Tensors" begin
        @test_throws DomainError Tensor([1, 2], VectorSpace(3), VectorSpace(2))
        @test_throws DomainError Tensor([1, 2], VectorSpace(3))
        @testset "DeltaTensor" begin
            @test DeltaTensor(V2, V2').shape == (V2, V2')
            @test_throws DomainError DeltaTensor(V2, VectorSpace(2))
            @test DeltaTensor().shape == ()
            @test_throws DomainError DeltaTensor(V2)
        end
    end
    @testset "Indices" begin
        @test_throws DomainError FixedIndex(VectorSpace(4), 5)
        @test_throws DomainError FixedIndex(VectorSpace(), 5)
    end
    @testset "Trigonometry" begin
        @test Base.sin(Tensor(1)) isa SineOperation
        @test Base.sin(Tensor(3)).children[1] == Tensor(3)
        @test Base.cos(VariableTensor("x")) isa CosineOperation
        @test Base.tan(B[1, 2]) isa TangentOperation
        @test_throws MethodError Base.sin(A)
        @test_throws MethodError Base.cos(A)
        @test_throws MethodError Base.tan(A)
        @test Base.asin(Tensor(.4)) isa ArcsineOperation
        @test_throws DomainError Base.asin(Tensor(4))
        @test Base.acos(Tensor(.4)) isa ArccosineOperation
        @test_throws DomainError Base.acos(Tensor(4))
        @test Base.atan(Tensor(40)) isa ArctangentOperation
    end
    @testset "Differentation" begin
        @test differentiate(a, a) == ConstantTensor(1)
        @test differentiate(a, VariableTensor("z")) == ZeroTensor()
        @test differentiate(a + ZeroTensor(), a) == ConstantTensor(1)
        @test differentiate(a + a, a) == ConstantTensor(1) + ConstantTensor(1)
        @test differentiate(a + 3*a, a) == ConstantTensor(1) + ConstantTensor(3)
        @test differentiate(1 / a, a) == (-1 / (a^2))
        @test differentiate(Base.sin(a), a) == Base.cos(a)
        @test differentiate(Base.cos(a), a) == -Base.sin(a)
        @test differentiate(Base.sin(a * VariableTensor("z")), a) == Base.cos(a * VariableTensor("z")) * VariableTensor("z")
        @test differentiate(Base.tan(Base.cos(a)), a) == - Base.sin(a) / (Base.cos(Base.cos(a))^2)
        @test differentiate(Base.asin(6*a), a) == 6 / Base.sqrt(Tensor(1) - (6a)^2)
        @test differentiate(Base.acos(6*a), a) == - Tensor(6) / Base.sqrt(Tensor(1) - (6a)^2)
        @test differentiate(Base.atan(6*a), a) == 6 / (Tensor(1) + (6a)^2)
        @test divergence(Tensor([b * sin(a), a * sin(b)], V2), a, b) == b * cos(a) + a * cos(b)
        @test divergence(Tensor([a, a], V2) + Tensor([b, b], V2), a, b) == Tensor(1) + Tensor(1)
        @test divergence(Tensor([43, a * sin(b)], V2), a, b) == a * cos(b)
        @test divergence(a * b * Tensor([5, 1], V2), a, b) == 5b + a
        @test divergence(a * (Tensor([0, 5b], V2) + Tensor([2b, 0], V2)), a, b) == 5a + 2b
        @test divergence(a * ConstantTensor(1, V2), a, b) == ConstantTensor(1)
        @test_throws DimensionMismatch divergence(Tensor([1, 2, 3], V3), a, b)
    end
    @testset "TreeVisitor" begin
        @testset "Update children" begin
            @test updatechildren(A + B, A, B, B).children == (A, B, B)
            @test updatechildren(A⊗B, A, E).children == (A, E)
            @test updatechildren(A[x, 1], A, fixedindices).children == (A, fixedindices,)
            @test updatechildren(A[x, y]*D[y', z], (B[x, y]*D[y', z]).children...) == B[x, y]*D[y', z]
            @test updatechildren(sin(VariableTensor("z")), Tensor(2)) == sin(Tensor(2))
            @test updatechildren(cos(VariableTensor("z")), Tensor(2)) == cos(Tensor(2))
            @test updatechildren(tan(VariableTensor("z")), Tensor(2)) == tan(Tensor(2))
        end
    end
    @testset "Assignment" begin
        @testset "Vector spaces" begin
            @test assign(D, Vi=>RnSpace(2)) == VariableTensor("D", V2', RnSpace(2))
            @test assign(ZeroTensor(Vi'), Vi'=>RnSpace(2)) == ZeroTensor(RnSpace(2))
            @test assign(ConstantTensor(a, Vi'), Vi'=>RnSpace(2)) == ConstantTensor(a, RnSpace(2))
            @test assign(DeltaTensor(Vi', Vi), Vi=>V2) == DeltaTensor(V2', V2)
            @test_throws DomainError assign(D, V3=>Vi)
            @test assign(D[1, z], Vi=>RnSpace(2)).freeindices == (FreeIndex(RnSpace(2), "z"),)
            @test assign(componenttensor(D[y', z], z, y'), Vi=>RnSpace(2)).shape == (RnSpace(2), V2')
        end
        @testset "Tensors" begin
            @test assign(A, A=>B) == B
            # TODO Create a better node equality so that strings are not needed
            @test string(assign((A⊗C)[2], A=>B)) == string((B⊗C)[2])
            @test_throws DomainError assign(D, D=>B)
            @test assign(A, A=>ConstantTensor(2, A.shape...)) == ConstantTensor(2, A.shape...)
            @test togem(DeltaTensor(V2, V2')).array == reshape([1, 0, 0, 1], (2, 2))
        end
        @testset "Variables" begin
            @test assign(a, a=>4) == ConstantTensor(4)
            @test assign(Tensor([a, 6], RnSpace(2)), a=>4).value == [4, 6]
            @test assign(Tensor([cos(a), a], V2), a => 0).value == [cos(ZeroTensor()), 0]
            @test assign(G, a=>0).value == [cos(ZeroTensor()), 4*sin(ZeroTensor()), 4, 7im]
        end
        @testset "Indices" begin
            @test assign(Indices(x), x=>1) == Indices(FixedIndex(V3, 1))
            @test string(assign(A[x], x=>1)) == string(A[1])
            @test string(assign(A[x, y], x=>1)) == string(A[1, y])
        end
    end
    @testset "Gem" begin
        @testset "Tensors" begin
            @test togem(B) == gem.Literal(fill(1.2, (3, 2)))
            @test togem(ConstantTensor(1, V3, V2)) == gem.Literal(fill(1, (3,2)))
            @test togem(H) == gem.Variable("H", (2, 5))
            @test togem(Z) == gem.Zero((3, 2))
        end
        @testset "Indices" begin
            @test togem(Indices(y, z, FixedIndex(V2, 2))) == togem(Indices(y, z, FixedIndex(V2, 2)))
            @test togem(H[y', 1]).free_indices == togem(Indices(y'))
            @test togem(H[y', 1]).children == (togem(H),)
            @test togem(H[y', 1]).multiindex == togem(Indices(y', FixedIndex(RnSpace(5), 1)))
            @test togem(Tensor(collect(1:5), RnSpace(5))[2]) == gem.Literal(2)
        end
        @testset "componenttensor" begin
            @test togem(A[x]).shape == (2,)
            @test togem(A[x]).free_indices == togem(Indices(x))
            @test togem(componenttensor(A[x],x)).shape == (2,3)
            @test togem(componenttensor(A[x],x)).free_indices == ()
        end
        @testset "Addition" begin
            @test togem(A + B).shape == (3, 2)
            @test isinst(togem(A + B), gem.ComponentTensor)
            @test isinst(togem(A + B).children[1], gem.Sum)
            @test togem(A + B + A).shape == (3, 2)
            @test togem(A + B).children[1].children[1].children[1] == togem(A)
            @test togem(A + B).children[1].children[2].children[1] == togem(B)
            @test togem(B + A + A).children[1].children[1].children[1].children[1] == togem(B)
            @test Set(togem(A[x, y] + B[1,2]).free_indices) == Set(togem(Indices(x, y)))
        end
        @testset "Product" begin
            @test togem(A⊗B).shape == (3, 2, 3, 2)
            @test togem(A⊗B).children[1].children[2].children[1] == togem(B)
            @test isinst(togem(A⊗B), gem.ComponentTensor)
            @test isinst(togem(A⊗B).children[1], gem.Product)
        end
        @testset "Division" begin
            @test isinst(togem(A[1,2]/3), gem.Division)
            @test togem(B[1,2] / 2) == gem.Literal(0.6)
            @test isinst(togem(A / VariableTensor("h")), gem.ComponentTensor)
            @test isinst(togem(A / VariableTensor("h")).children[1], gem.Division)
        end
        @testset "Index sum" begin
            @test togem(A[x, y]⊗F[y', 1]).free_indices == togem(Indices(x))
            @test togem(A[x, y]⊗F[y']).free_indices == togem(Indices(x))
            @test isinst(togem(A[x, y]⊗F[y']), gem.ComponentTensor)
            @test isinst(togem(A[x, y]⊗F[y']).children[1], gem.IndexSum)
            @test togem(A[x, y]⊗F[y']).free_indices == togem(Indices(x))
        end
        @testset "Trigonometry" begin
            @test togem(sin(a)).name == "sin"
            @test togem(cos(a)).name == "cos"
            @test togem(tan(a)).name == "tan"
            @test togem(asin(a)).name == "asin"
            @test togem(acos(a)).name == "acos"
            @test togem(atan(a)).name == "atan"
        end
    end
    # These tests only check that no errors are occurring, they do not check correctness.
    @testset "Create kernel" begin
        @test Kernel(A + B).knl !== nothing
        # This does not pass, but maybe it shouldn't?
        # @test Kernel(A[x, y]⊗F[y']) !== nothing
        @test Kernel((B + A + A)[1]).knl !== nothing
        @test Kernel(A[1, y]⊗F[y']).knl !== nothing
        @test Kernel(componenttensor((A + Z)[1, y]⊗F[y', x], x)).knl !== nothing
        @test Kernel(Tensor([a 1 ; 2 3], RnSpace(2), RnSpace(2))).knl !== nothing
    end
    @testset "execute" begin
        @test execute(B) == fill(1.2, (1, 3, 2))
        @test execute(B + 5*B) == fill(1.2 + 5*1.2, (1, 3, 2))
        @test execute(Kernel(A), Dict("A"=>fill(5.4, (3, 2)))) == fill(5.4, (1, 3, 2))
        @test execute(Kernel(A[1, y]⊗H[y']), Dict("A"=>fill(1.0, (3,2)), "H"=>fill(-1.0, (5,2)))) == fill(-2.0, (1,5))
        @test execute(I, "a"=>1.0) == reshape([0.0 + i for i in 1:12],(1,2,2,3))
        @test execute(Tensor([sin(a), cos(a), tan(a)], V3), "a"=>1.0) == reshape([sin(1), cos(1), tan(1)], (1, 3))
        @test execute(Tensor([asin(a), acos(a), atan(a)], V3), "a"=>0.5) == reshape([asin(0.5), acos(0.5), atan(0.5)], (1, 3))
        @test execute(A / a, "A"=>fill(3.0, (3, 2)), "a"=>1.5) == fill(2.0, (1, 3, 2))
        @test execute(trace(Tensor([1 2 ; 3 4], V2, V2'))) == [5.0]
        @test execute(sqrt(a), "a"=>4.0) == [2.0]
        @test execute(a^a, "a"=>3.0) == [27.0]
    end
    @testset "find variables" begin
        @test findvariables((A[x, y] + D[y', z])⊗B) == Set{VariableTensor}([A, D])
    end
end

max_stars_repo_name: Sohva/oolong

max_issues_repo_name: Sohva/oolong

max_forks_repo_name: oolongtensor/oolong


text:
using CSV
using MLJFair
using DataFrames

"
Utility function to read csv file of job data and convert columns to categorical
"
function jobdata()
    fpath = joinpath(@__DIR__, "..", "..", "data", "jobs.csv")
    data = DataFrame!(CSV.File(fpath))
    categorical!(data, names(data)[1:6])
    return data
end

"
Utility function for tests to create fairness Tensor for jobs.csv
The fairness tensor will be created on the basis of the column Job Type.
This column has 3 different values for job types
"
function job_fairtensor()
    data = jobdata()
    return fair_tensor(data[!, :Pred], data[!, :Class], data[!, names(data)[4]])
end


text:
import DandelionWebSockets: wsconnect, on_text, on_binary
export on_text, on_binary

struct GDAXWebSocketClient{T <: AbstractGDAXMessageHandler} <: DandelionWebSockets.WebSocketHandler
    websocket_feed::String
    client::WSClient
    subscription::Dict{String, Any}
    events_handler::T

    function GDAXWebSocketClient(user::GDAXUser,
                                 subscription::Dict{String, Any},
                                 events_handler::T) where {T <: AbstractGDAXMessageHandler}

        if isregistered(user)
            merge!(subscription, getHeaders(user, "GET", "/users/self/verify", Dict()))
        end

        return new{T}(user.ws_feed, WSClient(ponger = DandelionWebSockets.Ponger(3.0; misses=90)), subscription, events_handler)
    end
end

function wsconnect(this::GDAXUser, subscription::Dict{String, Any}, events_handler::T) where {T <: AbstractGDAXMessageHandler}
    client = GDAXWebSocketClient(this, subscription, events_handler)
    connect(client)
    subscribe(client)
    return client
end

function connect(this::GDAXWebSocketClient)::GDAXWebSocketClient
    flag = DandelionWebSockets.wsconnect(this.client, URI(this.websocket_feed), this.events_handler)
    if !flag
        ws = this.websocket_feed
        throw(ErrorException("Unable to connect to $ws"))
    end

    return this
end

function subscribe(this::GDAXWebSocketClient)
    DandelionWebSockets.send_text(this.client, JSON.json(this.subscription))
end

function unsubscribe(this::GDAXWebSocketClient)
    unsubscription = this.subscription
    unsubscription["type"] = "unsubscribe"
    DandelionWebSockets.send_text(this.client, JSON.json(unsubscription))
end

function onWSMessage(this::AbstractGDAXMessageHandler, msg::Dict{String, X}) where {X <: Any}
    T = typeof(this)
    throw(ErrorException("Method `onMessage` is not implemented by $T"))
end

function on_text(this::T, str::String)::Void where {T <: AbstractGDAXMessageHandler}
    msg = try
        JSON.parse(str)
    catch exc
        Dict("type" => "parse_error",
             "text" => str,
             "exception" => exc)
    end

    onWSMessage(this, msg)
    return nothing
end

function on_binary(this::T, msg::Vector{UInt8}) where {T <: AbstractGDAXMessageHandler}
    @printf("[GDAXClient] Incoming binary data (this is unexpected)")
end

using .SlowOrderBook

struct GDAXOrderBooks <: AbstractGDAXMessageHandler
    books::Dict{String, OrderBook}
    condition::Condition
    function GDAXOrderBooks(products::Vector{String}, condition::Condition)
        b = Dict{String, OrderBook}()
        [b[x] = OrderBook() for x in products]
        return new(b, condition)
    end
end

import Base: getindex
export getindex

function getindex(this::GDAXOrderBooks, product_id::String)
    return this.books[product_id]
end

function onWSMessage(this::GDAXOrderBooks, msg::Dict{String, X}) where {X <: Any}
    if msg["type"] == "snapshot"
        book = this.books[msg["product_id"]]
        for asks in msg["asks"]
            price = parse(Float64, asks[1])
            lots = parse(Float64, asks[2])
            book["sell"][price] = lots
        end

        for bids in msg["bids"]
            price = parse(Float64, bids[1])
            lots = parse(Float64, bids[2])
            book["buy"][price] = lots
        end

        book.now.data = now(Dates.UTC)
    elseif msg["type"] == "l2update"
        book = this.books[msg["product_id"]]
        for item in msg["changes"]
            side = item[1]
            price = parse(Float64, item[2])
            lots = parse(Float64, item[3])

            book[side][price] = lots
        end

        book.now.data = DateTime(msg["time"], "yyyy-mm-ddTHH:MM:SS.sss\\Z")
        notify(this.condition)
    end
end

max_stars_repo_name: DatName/gdax-julia

max_issues_repo_name: DatName/gdax-julia

max_forks_repo_name: DatName/gdax-julia