text:
import DataFrames: DataFrame

"""

`type_array(arr)`

`type_array` finds the most constraining type of the elements of an array `arr`, and converts the array element type.
Sometimes, an array is given type `Array{Any}`, even though the components are all Float64, for example.
`type_array` will convert the type into the most constraining one.

##### Arguments
* `arr::AbstractArray`: an abstract array whose element type will be constrained by `type_array`.

##### Returns
An array with the same elements as in `arr`, but the element type has been constrained just enough to contain all elements.

##### Examples

```julia
julia> type_array(Any[1, 3.0, 2])
3-element Array{Float64,1}:
 1.0
 3.0
 2.0

julia> type_array(Any[1, 3.0, 'x'])
3-element Array{Any,1}:
 1
 3.0
  'x'
```

"""
function type_array end

type_array(arr::FloatNAArray) = arr

type_array{F,N,A<:FloatNAArray}(x::AbstractArrayWrapper{F,N,A}) = x

type_array{T,N}(x::AbstractArray{Nullable{T},N}) = x

type_array{T}(x::AbstractArray{T}) = begin
  if isempty(x) # || (isa(x, AbstractArrayWrapper) && isa(x.a, FloatNAArray))
    x
  else
    #eltypex = T #eltype(x)
    # return x intact if it is already an array of a concrete type.
    # but Nullable{T} is also a concrete type... Let's drop the concrete type check part in this case.
    # well, not having this makes the code more than twice faster...
    # if !isa(T,Union) && !(T<:Nullable) && isempty(subtypes(T))
    #  x
    if !isdefined(x, 1)
      x
    else
      common_type = typeof(x[1])
      # if the first element is already the same as eltype(x), there is no chance of further restriction.
      # still, don't understand why == does not work.
      if common_type === T
        return x
      end
      type_changed = false
      for elem in x
        elemtype = typeof(elem)
        if !type_changed && elemtype !== T
          type_changed = true
        end
        if !isa(elem, common_type)
          common_type = promote_type(common_type, elemtype)
        end
      end
      if type_changed
        result = similar(x, common_type)
        copy!(result, x)
        result
      else
        x
      end
    end
  end
end
type_array(arr::DictArray) = arr
type_array(arr::LabeledArray) = arr
type_array(arr::DataFrame) = arr
type_array(arr::EnumerationArray) = arr

# used in join functions. broadcast an array of dimensions N1x...xNn into
# size(front_dims) x N1 x ... x Nn x size(back_dims).
expand_dims(arr::LabeledArray, front_dims::NTuple{TypeVar(:M),AbstractArray}, back_dims::NTuple{TypeVar(:N),AbstractArray}) = begin
  front_size = map(length,front_dims)
  back_size = map(length,back_dims)
  data = expand_dims(arr.data, front_size, back_size)
  axes = (front_dims..., arr.axes..., back_dims...)
  LabeledArray(data, axes)
end

expand_dims(arr::LabeledArray, front_dims::Tuple{}, back_dims::Tuple{}) = arr

expand_dims(arr::DictArray, front_size::NTuple{TypeVar(:U),Int}, back_size::NTuple{TypeVar(:V),Int}) = begin
  DictArray(LDict(arr.data.keys, map(x->expand_dims(x, front_size, back_size), arr.data.values)))
end

# There are M nested for loops generated by @generate.
# I need to pass the variable M as a number, so I define it as a mutable function.
@generated expand_dims{N,M,L}(arr::AbstractArray{TypeVar(:T),N}, front_size::NTuple{M,Int}, back_size::NTuple{L,Int}) = begin
  tgt_ndims = N + M + L
  srcindices = ntuple(n->symbol("i_",n+M), N)
  srcexp = Expr(:ref, :arr, srcindices...)
  quote
    data = similar(arr, (front_size...,size(arr)...,back_size...))
    @nloops $tgt_ndims i data begin
      @nref($tgt_ndims,data,i) = $(srcexp)
    end
    data
  end
end

# used internally as a default function to create a new field name.
generic_fieldname(i::Int) = symbol('x', i)
is_generic_fieldname(fieldname) = isa(fieldname, Symbol) && startswith(string(fieldname), "x") && try
  parse(Int, string(fieldname)[2:end])
  true
catch err
  if isa(err, ArgumentError)
    false
  else
    rethrow()
  end
end

# Create a new field name not currently used as a field name in the LabeledArray.
# Creates only a symbol.
create_additional_fieldname(arr::LabeledArray, tracker=nothing, create_fieldname::Function=generic_fieldname) = begin
  axesnames = map(filter(x->isa(x,DictArray), arr.axes)) do axis
    axis.data.keys
  end
  allnames = if isa(arr.data, DictArray)
    [arr.data.data.keys;axesnames...]
  else
    [axesnames...;]
  end
  i = 1
  while(true)
    label = create_fieldname(i)
    if !in(label, allnames)
      if tracker == nothing
        return label
      elseif !in(label, tracker)
        push!(tracker, label)
        return label
      end
    end
    i += 1
  end
end

"""

`collapse_axes(arr::AbstractArray, front_dim::Integer, back_end::Integer)`

Collapse front_dim to back_dim dimensions into one.

##### Arguments

* `arr` : an array
* `front_dims` : the starting direction to collapse.
* `back_dims` : the end direction to collapse.

The result is an array whose elements along front_dims to back_dims are all flattened into one dimension.
If `arr` is a LabeledArray, all the labels along the flattened direction are combined together.

##### Examples

```julia
julia> collapse_axes(darr(a=reshape(1:40, 2,4,5), b=reshape(11:50, 2,4,5)), 1, 2)
8 x 5 DictArray

a b  |a  b  |a  b  |a  b  |a  b
-----+------+------+------+------
1 11 |9  19 |17 27 |25 35 |33 43
2 12 |10 20 |18 28 |26 36 |34 44
3 13 |11 21 |19 29 |27 37 |35 45
4 14 |12 22 |20 30 |28 38 |36 46
5 15 |13 23 |21 31 |29 39 |37 47
6 16 |14 24 |22 32 |30 40 |38 48
7 17 |15 25 |23 33 |31 41 |39 49
8 18 |16 26 |24 34 |32 42 |40 50
```

"""
function collapse_axes end

collapse_axes(arr::AbstractArray, front_dim::Integer, back_dim::Integer) = begin
  sizearr = collect(size(arr))
  newsize = [sizearr[1:front_dim-1];prod(sizearr[front_dim:back_dim]);sizearr[back_dim+1:end]]
  reshape(arr, (newsize...))
end

# No need to do song and dance for 1 dimensinoal arrays.
collapse_axes(arr::AbstractArray{TypeVar(:V),1}, front_dim::Integer, back_dim::Integer) = begin
  @assert(front_dim == 1)
  @assert(back_dim == 1)
  arr
end
collapse_axes(arr::DictArray{TypeVar(:K),1}, front_dim::Integer, back_dim::Integer) = begin
  @assert(front_dim == 1)
  @assert(back_dim == 1)
  arr
end
collapse_axes(arr::LabeledArray{TypeVar(:T),1}, front_dim::Integer, back_dim::Integer) = begin
  @assert(front_dim == 1)
  @assert(back_dim == 1)
  arr
end

collapse_axes(arr::DictArray, front_dim::Integer, back_dim::Integer) =
  DictArray(LDict(arr.data.keys, map(v->collapse_axes(v, front_dim, back_dim), arr.data.values)))

collapse_axes(arr::LabeledArray, front_dim::Integer, back_dim::Integer) = begin
  ndimsarr = ndims(arr)
  sizearr = size(arr)
  newdata = collapse_axes(arr.data, front_dim, back_dim)
  arraxes = collect(arr.axes)
  newaxes = cell(ndimsarr - (back_dim-front_dim))
  newaxes[1:front_dim-1] = arraxes[1:front_dim-1]
  newaxes[front_dim+1:end] = arraxes[back_dim+1:end]
  collapsed_axes_vec = map(filter(d->!isa(arraxes[d],DefaultAxis), front_dim:back_dim)) do d
    inner_dims = prod(sizearr[front_dim:d-1])
    outer_dims = prod(sizearr[d+1:back_dim])
    repeat(arraxes[d], inner=collect(inner_dims), outer=collect(outer_dims))
  end
  collapsed_axes = if isempty(collapsed_axes_vec)
    DefaultAxis(prod(map(length, arraxes[front_dim:back_dim])))
  else
    tracker = []
    DictArray(merge(map(x->isa(x,DictArray) ? x.data : LDict(create_additional_fieldname(arr,tracker)=>x),
                           collapsed_axes_vec)...))
  end
  newaxes[front_dim] = collapsed_axes
  LabeledArray(newdata, (newaxes...))
end

collapse_axes(arr::AbstractArray) = collapse_axes(arr, 1, ndims(arr))

"""

`@rap(args...)`

Apply right-to-left evaluation order to the arguments.
An argument having an underscore symbol except the last one is translated into the function `x`->(that expression replacing _ by `x`).

##### Examples

```julia
julia> @rap _+3 5
8

julia> @rap _*2 x->x+1 10
22

julia> @rap (_ .* 2) reverse @nalift [1,2,NA,4,5]
5-element DataCubes.AbstractArrayWrapper{Nullable{Int64},1,Array{Nullable{Int64},1}}:
 Nullable(10)
 Nullable(8)
 Nullable{Int64}()
 Nullable(4)
 Nullable(2)
```
"""
macro rap(args...)
  esc(recursive_rap(args...))
end

recursive_rap(expr)= expr
recursive_rap(expr0, exprs...) = begin
  if does_underscore_exist(expr0)
    if :head in fieldnames(expr0) && expr0.head == Symbol("'")
      replaced_expr = replace_underscore(expr0, :x)
      quote
        mapna(x->$replaced_expr, $(recursive_rap(exprs...)))
      end
    else
      replace_underscore(expr0, recursive_rap(exprs...))
    end
  else
    if :head in fieldnames(expr0) && expr0.head == Symbol("'")
      quote
        mapna($expr0, ($(recursive_rap(exprs...))))
      end
    else
      quote
        $expr0($(recursive_rap(exprs...)))
      end
    end
  end
end
quote_symbol(s) = if isa(s, Symbol); QuoteNode(s) else s end

does_underscore_exist(expr) = if expr == :_
  true
elseif :head in fieldnames(expr)
  does_underscore_exist(expr.head) || any(Bool[does_underscore_exist(arg) for arg in expr.args ])
else
  false
end

replace_underscore(expr, repby) = if expr == :_
  repby
elseif :head in fieldnames(expr)
  Expr(replace_underscore(expr.head, repby), map(x->replace_underscore(x, repby), expr.args)...)
else
  expr
end

is_kw(arg) = :head in fieldnames(arg) && (arg.head==:kw || arg.head==:(=) || arg.head==:(=>))

read_kws(args) = begin
  map(args) do arg
    if !is_kw(arg)
      throw(ArgumentError("not a key=value type argument: $arg"))
    end
    key = if arg.head == :kw || arg.head == :(=); quote_symbol(arg.args[1]) else arg.args[1] end
    value = arg.args[2]
    quote
      $key => @nalift($(esc(value)))
    end
  end
end

"""

##### Description
Replace axes with another fields.
The args are a list of pairs of the form (integer for the axis index) => new axes fields for this axis.
Only the first elements (`arr[:,...,:,1,:,...,:]`) will be taken. That is, if the underlying data array is 2 dimensional, and
you want to use the field `column1` as a new key for the 1st axis, `column1[:,1]` will be used as the new axis.
e.g. `replace_axes(labeled_array, 1=>[:c1,:c2], 3=>[:c3])`

"""
function replace_axes end

replace_axes(arr::LabeledArray) = arr

"""
#### `replace_axes(arr::LabeledArray, args...)`

##### Arguments
* `arr`: an input `LabeledArray`. Note the data part of `arr` and the axis to replace need to be `DictArray`s.
* `args...`: arguments in `args` are of the form `i=>[f1,f2,...]` for some integer `i` and field names `f*`.

##### Returns
For each `i=>[f1,f2,...]` in `args`, the `i`th axis is replaced by the fields `f1`, `f2`, ....
A nonarray value in `i=>f1` will be automatically promoted to an array with one element: `[f1]`.
If the key `i` is omitted, the smallest relevant axis number will be automatically assigned.
Only the first elements will be taken. For example, if the underlying data array is 2 dimensional, and if you want to use some field for the 1st axis, `[:,1]` components will be used.
The original axis becomes the data part of `LabeledArray` after properly broadcast.
If the field name array is null for an argument in `args` (`i=>[]`), the corresponding axis will be `DefaultAxis`.

##### Examples

```julia
julia> t = larr(a=[1 2 3;4 5 6], b=['a' 'b' 'c';'d' 'e' 'f'], axis1=darr(k=[:x,:y]), axis2=["A","B","C"])
2 x 3 LabeledArray

  |A   |B   |C
--+----+----+----
k |a b |a b |a b
--+----+----+----
x |1 a |2 b |3 c
y |4 d |5 e |6 f

julia> replace_axes(t, 1=>[:a])      # == replace_axes(t, [:a]) == replace_axes(t, :a)
2 x 3 LabeledArray

  |A   |B   |C
--+----+----+----
a |b k |b k |b k
--+----+----+----
1 |a x |b x |c x
4 |d y |e y |f y

julia> replace_axes(t, 1=>[:a, :b])
2 x 3 LabeledArray

    |A |B |C
----+--+--+--
a b |k |k |k
----+--+--+--
1 a |x |x |x
4 d |y |y |y
```

"""
replace_axes(arr::LabeledArray, args...; index_counter::Int=1) = begin
  new_index_counter, axis_pair = replace_axes_inner_create_axis_pair(index_counter, args[1])
  replace_axes(replace_axes(arr, axis_pair; index_counter=1), args[2:end]...; index_counter=new_index_counter)
end

replace_axes(arr::LabeledArray, arg; index_counter::Int=1) = begin
  res = arr
  numdims = ndims(arr)
  tracker = []
  (_, (axis_index, axis_keys_raw)) = replace_axes_inner_create_axis_pair(index_counter, arg)
  axis_keys = replace_axes_inner_collect_if_necessary(axis_keys_raw)
  resaxis = res.axes[axis_index]
  axis = if isa(resaxis, DictArray)
    resaxis.data
  elseif isa(resaxis, DefaultAxis)
    LDict()
  else
    newkey = create_additional_fieldname(arr, tracker)
    axiswithkey = create_ldict_nocheck(newkey => resaxis)
    res = larr(res, symbol(:axis,axis_index) => create_dictarray_nocheck(axiswithkey))
    axiswithkey
  end
  orig_axis_keys = axis.keys
  new_axis = if isempty(axis_keys)
    DefaultAxis(size(arr)[axis_index])
  else
    DictArray((map(axis_keys) do axis_key
      if in(axis_key, orig_axis_keys)
        axis_key => axis[axis_key]
      else
        coords = ntuple(numdims) do i
          if i == axis_index
            Colon()
          else
            1
          end
        end
        axis_key => wrap_array(collect(selectfield(res, axis_key)[coords...]))
      end
    end)...)
  end
  resdict = LDict{Any,Any}()
  for (k,v) in res.data.data
    if !in(k, axis_keys)
      resdict[k] = v
    end
  end
  for (k,v) in axis
    if !in(k, axis_keys)
      resdict[k] = selectfield(res, k)
    end
  end
  axes = ntuple(numdims) do d
    if d == axis_index
      new_axis
    else
      res.axes[d]
    end
  end
  LabeledArray(DictArray(resdict...), axes)
end

replace_axes_inner_collect_if_necessary(axis_keys_raw) = [axis_keys_raw]
replace_axes_inner_collect_if_necessary(axis_keys_raw::AbstractArray) = axis_keys_raw

replace_axes_inner_create_axis_pair(index_counter::Int, arg::Pair) = (arg[1]+1, arg)
replace_axes_inner_create_axis_pair(index_counter::Int, arg::NTuple{2}) = (arg[1]+1, arg)
replace_axes_inner_create_axis_pair(index_counter::Int, arg) = (index_counter+1, index_counter => arg)

"""

`flds2axis(arr::LabeledArray [; axisname=nothing, fieldname=nothing])`

Create another dimension using the field values of the data of a LabeledArray.

##### Arguments

* `arr` : a `LabeledArray`.
* `axisname` (optional) : the name of the new axis.
* `fieldname` (optional) : the name of the new field name. If not specified, the resulting `LabeledArray` will have a normal `AbstractArray` and not a `DictArray` as its data.

##### Returns

A new `LabeledArray` which has one higher dimensions than the input `arr`.
The field names become the elements of the new last axis, after wrapped by `Nullable`.
If `axisname` is provided, the new axis becomes a `DictArray` with that field name.
Otherwise, the new axis will be a normal array.
If `fieldname` is provided, the new data of the return `LabeledArray` is a `DictArray` with that field name.
Otherwise, the new data will be a normal array.

##### Examples

```julia
julia> t = larr(a=[1,2,3], b=[:x,:y,:z])
3 LabeledArray

  |a b
--+----
1 |1 x
2 |2 y
3 |3 z

julia> flds2axis(t, axisname=:newaxis, fieldname=:newfield)
3 x 2 LabeledArray

newaxis |a        |b
--------+---------+---------
        |newfield |newfield
--------+---------+---------
1       |1        |x
2       |2        |y
3       |3        |z
```

"""
function flds2axis end

flds2axis(arr::LabeledArray; axisname=nothing, fieldname=nothing) = begin
  ldict = arr.data.data
  eltypes = [eltype(v) for (_,v) in ldict]
  common_eltype = promote_type(eltypes...)
  resultdata = wrap_array(similar(arr.data, common_eltype, (size(arr)...,length(ldict))))
  colons = ntuple(_ -> Colon(), ndims(arr))
  for i in 1:length(ldict)
    resultdata[colons..., i] = ldict.values[i]
  end
  newaxis = if axisname == nothing
    nalift(ldict.keys)
  else
    DictArray(axisname => nalift(ldict.keys))
  end
  newdata = if fieldname == nothing
    resultdata
  else
    DictArray(fieldname => resultdata)
  end
  LabeledArray(newdata, (arr.axes..., newaxis))
end

"""

`axis2flds(arr::LabeledArray (; name_collapse_function=..., default_axis_value=nothing)`

Collapse a dimension of a [`LabeledArray`](#type__labeledarray.1), making the axis along that direction as field names.

##### Arguments

* `larr` : a LabeledArray.
* `name_collapse_function` (optional keyword) : a function to combine the axis label and the column name. By default, it concatenates the names with '_' inserted in between.
* `default_axis_value` (optional keyword) : a default value to be used when an axis element is null. If `nothing` (by default), an exception will raise.

##### Examples

```julia
julia> t = larr(reshape(1:10,5,2), axis1=darr(k=['a','b','c','d','e']), axis2=darr(r1=[:M,:N],r2=["A","A"]))
5 x 2 LabeledArray

r1 |M |N
r2 |A |A
---+--+---
k  |  |
---+--+---
a  |1 |6
b  |2 |7
c  |3 |8
d  |4 |9
e  |5 |10

julia> axis2flds(t)
5 LabeledArray

k |M_A N_A
--+--------
a |1   6
b |2   7
c |3   8
d |4   9
e |5   10

julia> axis2flds(t, name_collapse_function=x->join(x, "*"))
5 LabeledArray

k |M*A N*A
--+--------
a |1   6
b |2   7
c |3   8
d |4   9
e |5   10

julia> m = @larr(reshape(1:10,5,2), axis1[k=['a','b','c','d','e']], axis2[:M,NA])
5 x 2 LabeledArray

  |M |
--+--+---
k |  |
--+--+---
a |1 |6
b |2 |7
c |3 |8
d |4 |9
e |5 |10

julia> axis2flds(m, default_axis_value="N/A")
5 LabeledArray

k |M N/A
--+------
a |1 6
b |2 7
c |3 8
d |4 9
e |5 10
```

"""
function axis2flds end

axis2flds(larr::LabeledArray; name_collapse_function=default_fieldname_collapse_function, default_axis_value=nothing) = begin
  lastaxis = larr.axes[ndims(larr)]
  lastaxiskeys = if isa(lastaxis, DictArray)
    # null value is not allowed for an axis value if you want to collapse it.
    map(lastaxis.data.values...) do values...
      (map(x->get_nullable_value(x, default_axis_value), values)...)
    end
  else
    map(x->(get_nullable_value(x, default_axis_value),), lastaxis)
  end
  colons = ntuple(_ -> Colon(), ndims(larr) - 1)
  newdictvec = if isa(larr.data, DictArray)
    [LDict(map(x->name_collapse_function(lastaxiskeys[i], x), larr.data.data.keys), larr.data[colons..., i].data.values) for i in 1:last(size(larr))]
  else
    [LDict(name_collapse_function(lastaxiskeys[i]) => larr.data[colons..., i]) for i in 1:size(larr,ndims(larr))]
  end
  newdata = darr(merge(newdictvec...)...)
  newaxes = larr.axes[1:ndims(larr)-1]
  LabeledArray(newdata, newaxes)
end

get_nullable_value(x::Nullable, ::Void) = if x.isnull
  throw(ArgumentError("an element in the last axis is null. Consider using an optional argument default_axis_value=..."))
else
  x.value
end
get_nullable_value(x::Nullable, default_axis_value) = x.isnull ? default_axis_value : x.value

const default_name_delimeter = '_'
default_fieldname_collapse_function(axistuple, field) = begin
  symbol(map(x->string(x, default_name_delimeter), axistuple)..., string(field))
end
default_fieldname_collapse_function(axistuple) = if length(axistuple) == 1
  axistuple[1]
else
  symbol(join(axistuple, default_name_delimeter))
end

"""

`mapvalues(f::Function, xs...)`

Apply a function `f` to each tuple constructed from `xs`. Each element of `xs` will be of type `LDict`/`DictArray`/`LabeledArray`. The input tuple is constructed by combining the element in the same position in each `xs` (or the constant value in case it is not an array).

##### Returns

For each element `x` in `xs`,

* If `x` is `LDict`, `f` is applied to each value and the result is again `LDict` with the same keys and the new values.
* If `x` is `DictArray`, `f` is applied to each field. The return value will be `DictArray` if the return value of `f` is also an `AbstractArray`. Otherwise, an `LDict` will be returned.
* If `x` is `LabeledArray`, `mapvalues(f, _)` is applied to the base of `x`. The return value will be `LabeledArray` with the same axes if the return value of `f` is also an `AbstractArray`. Otherwise, an `LDict` will be returned.

##### Examples

```julia
julia> mapvalues(x->x+1, LDict(:a=>1, :b=>2))
DataCubes.LDict{Symbol,Int64} with 2 entries:
  :a => 2
  :b => 3

julia> mapvalues(x->x .+ 1, darr(a=[1,2,3], b=[4,5,6]))
3 DictArray

a b
----
2 5
3 6
4 7

julia> mapvalues(x->x .+ 1, larr(a=[1,2,3], b=[4,5,6], axis1=[:m,:n,:p]))
3 LabeledArray

  |a b
--+----
m |2 5
n |3 6
p |4 7

julia> mapvalues(sum, darr(a=[1,2,3], b=[4,5,6]))
DataCubes.LDict{Symbol,Nullable{Int64}} with 2 entries:
  :a => Nullable(6)
  :b => Nullable(15)

julia> mapvalues(sum, larr(a=[1,2,3], b=[4,5,6], axis=[:m,:n,:p]))
DataCubes.LDict{Symbol,Nullable{Int64}} with 2 entries:
  :a => Nullable(6)
  :b => Nullable(15)

julia> mapvalues((x,y)->2x.*y,darr(a=[1 2 3]),darr(a=[4 5 6]))
1 x 3 DictArray

a |a  |a
--+---+---
8 |20 |36
```

"""
function mapvalues end

mapvalues(f::Function, arr::AbstractArray) = map(f, arr)
mapvalues(f::Function, arr::DictArray) = begin
  results = mapvalues(f, arr.data)
  if mapvalues_helper_allarrays_checker(results.values...)
    create_dictarray_nocheck(results)
  else
    results
  end
end
mapvalues(f::Function, arr::LabeledArray) = begin
  if isa(arr.data, DictArray)
    results = mapvalues(f, arr.data.data)
    if mapvalues_helper_allarrays_checker(results.values...)
      LabeledArray(create_dictarray_nocheck(results), arr.axes)
    else
      results
    end
  else
    LabeledArray(mapvalues(f, arr.data), arr.axes)
  end
end
mapvalues(f::Function, xs...) = begin
  allkeys = union(filter(x->!isempty(x), map(mapvalues_getkeys_helper, xs))...)
  reordered_xs = map(x->mapvalues_reorder_helper(allkeys,x), xs)
  isanydarr = any(x->isa(x,DictArray) || (isa(x,LabeledArray) && isa(peel(x),DictArray)), xs)
  labels = Nullable()
  for x in xs
    if isa(x, LabeledArray)
      if labels.isnull
        labels = Nullable(pickaxis(x))
      elseif labels.value != pickaxis(x)
        throw(ArgumentError("axes do not match."))
      end
    end
  end
  if isempty(allkeys)
    result = f(reordered_xs...)
    if labels.isnull
      result
    else
      LabeledArray(result, labels.value)
    end
  else
    result = map(reordered_xs...) do elems...
      f(elems...)
    end
    if mapvalues_helper_allarrays_checker(result...)
      if isanydarr
        if labels.isnull
          DictArray(allkeys, result)
        else
          LabeledArray(DictArray(allkeys, result), labels.value)
        end
      end
    else
      LDict(allkeys, result)
    end
  end
end

mapvalues_getkeys_helper(xs) = []
mapvalues_getkeys_helper(xs::LDict) = keys(xs)
mapvalues_getkeys_helper(xs::DictArray) = mapvalues_getkeys_helper(peel(xs))
mapvalues_getkeys_helper(xs::LabeledArray) = mapvalues_getkeys_helper(peel(xs))
mapvalues_reorder_helper(allkeys, xs) = isempty(allkeys) ? xs : fill(xs, length(allkeys))
mapvalues_reorder_helper(allkeys, xs::LDict) = keys(xs)==allkeys ? values(xs) : values(extract(xs, allkeys))
mapvalues_reorder_helper(allkeys, xs::DictArray) = mapvalues_reorder_helper(allkeys, peel(xs))
mapvalues_reorder_helper(allkeys, xs::LabeledArray) = mapvalues_reorder_helper(allkeys, peel(xs))

mapvalues_helper_allarrays_checker(::AbstractArray...) = true
mapvalues_helper_allarrays_checker(::Any...) = false

# permutedims only if necessary. Otherwise, return the argument itself.
# used mainly to reshuffle an array as an intermediate result, and
# when a copy is not necessarily required.
permutedims_if_necessary(arr::AbstractArray, perms) = ntuple(identity,ndims(arr))==perms ? arr : permutedims(arr, perms)
ipermutedims_if_necessary(arr::AbstractArray, perms) = ntuple(identity,ndims(arr))==perms ? arr : ipermutedims(arr, perms)

"""

Peel off a variable to see its underlying data.

* `peel(arr::DictArray)`: returns an `LDict` consisting of field name => field values array pairs.

* `peel(arr::LabeledArray)`: returns the underlying data, which can be a `DictArray` but can also be any `AbstractArray`.

* `peel(arr::EnumerationArray)`: returns the underlying index integers.

##### Examples

```julia
julia> peel(darr(a=[1,2,3], b=[:m,:n,:p]))
DataCubes.LDict{Symbol,DataCubes.AbstractArrayWrapper{T,1,A<:AbstractArray{T,N}}} with 2 entries:
  :a => [Nullable(1),Nullable(2),Nullable(3)]
  :b => [Nullable(:m),Nullable(:n),Nullable(:p)]

julia> peel(larr(a=[1,2,3], b=[:m,:n,:p], axis1=["X","Y","Z"]))
3 DictArray

a b
----
1 m
2 n
3 p

julia> peel(@enumeration([NA :x :y;:x :z NA]))
2x3 DataCubes.AbstractArrayWrapper{Int64,2,Array{Int64,2}}:
 0  1  3
 1  2  0

```

"""
function peel end

peel(arr::DictArray) = arr.data
peel(arr::LabeledArray) = arr.data
peel(arr::EnumerationArray) = arr.elems

"""

Pick fields from a `DictArray` or a `LabeledArray`.

* `pick(arr::DictArray, fieldname)`: returns the field value array corresponding to `fieldname`.
* `pick(arr::DictArray, fieldnames::AbstractArray)`: returns a `DictArray` whose field names are `fieldnames`.
* `pick(arr::DictArray, fieldnames::Tuple)`: returns a vector of field value arrays corresponding to `fieldnames`.
* `pick(arr::DictArray, fieldnames::...)` if there are more than 1 field name in `fieldnames`: returns a vector of field value arrays corresponding to the `fieldnames`.
* `pick(arr::LabeledArray, fieldname)`: returns the field value array corresponding to `fieldname`. If `fieldname` corresponds to a field in an axis, the field value array is broadcast appropriately.
* `pick(arr::LabeledArray, fieldnames::AbstractArray)`: returns a `DictArray` whose field names are `fieldnames`.
* `pick(arr::LabeledArray, fieldnames::Tuple)`: returns a vector of field value arrays corresponding to `fieldnames`.
* `pick(arr::LabeledArray, fieldnames::...)` if there are more than 1 field name in `fieldnames`: returns a vector of field value arrays corresponding to the `fieldnames`.

##### Examples

```julia
julia> pick(darr(a=[1,2,3], b=[:m,:n,:p]), :a)
3-element DataCubes.AbstractArrayWrapper{Nullable{Int64},1,Array{Nullable{Int64},1}}:
 Nullable(1)
 Nullable(2)
 Nullable(3)

julia> pick(darr(a=[1,2,3], b=[:m,:n,:p]), (:a,))
1-element Array{DataCubes.AbstractArrayWrapper{T,1,A<:AbstractArray{T,N}},1}:
 [Nullable(1),Nullable(2),Nullable(3)]

julia> pick(darr(a=[1,2,3], b=[:m,:n,:p]), :a, :b)
2-element Array{DataCubes.AbstractArrayWrapper{T,1,A<:AbstractArray{T,N}},1}:
 [Nullable(1),Nullable(2),Nullable(3)]
 [Nullable(:m),Nullable(:n),Nullable(:p)]

julia> pick(darr(a=[1,2,3], b=[:m,:n,:p]), (:a, :b))
2-element Array{DataCubes.AbstractArrayWrapper{T,1,A<:AbstractArray{T,N}},1}:
 [Nullable(1),Nullable(2),Nullable(3)]
 [Nullable(:m),Nullable(:n),Nullable(:p)]

julia> t = larr(a=[1 2;3 4;5 6], b=[:x :y;:z :u;:v :w], axis1=darr(k=["X","Y","Z"]), axis2=[:A,:B])
3 x 2 LabeledArray

  |A   |B
--+----+----
k |a b |a b
--+----+----
X |1 x |2 y
Y |3 z |4 u
Z |5 v |6 w

julia> pick(t, :a)
pic3x2 DataCubes.AbstractArrayWrapper{Nullable{Int64},2,Array{Nullable{Int64},2}}:
 Nullable(1)  Nullable(2)
 Nullable(3)  Nullable(4)
 Nullable(5)  Nullable(6)

julia> pick(t, :a, :k)
2-element Array{DataCubes.AbstractArrayWrapper{T,N,A<:AbstractArray{T,N}},1}:
 3x2 DataCubes.AbstractArrayWrapper{Nullable{Int64},2,Array{Nullable{Int64},2}}:
 Nullable(1)  Nullable(2)
 Nullable(3)  Nullable(4)
 Nullable(5)  Nullable(6)
 3x2 DataCubes.AbstractArrayWrapper{Nullable{ASCIIString},2,DataCubes.BroadcastAxis{Nullable{ASCIIString},2,DataCubes.AbstractArrayWrapper{Nullable{ASCIIString},1,Array{Nullable{ASCIIString},1}},DataCubes.DictArray{Symbol,2,DataCubes.AbstractArrayWrapper{T,2,A<:AbstractArray{T,N}},Nullable{T}}}}:
 Nullable("X")  Nullable("X")
 Nullable("Y")  Nullable("Y")
 Nullable("Z")  Nullable("Z")

julia> pick(t, (:a, :k))
2-element Array{DataCubes.AbstractArrayWrapper{T,N,A<:AbstractArray{T,N}},1}:
 3x2 DataCubes.AbstractArrayWrapper{Nullable{Int64},2,Array{Nullable{Int64},2}}:
 Nullable(1)  Nullable(2)
 Nullable(3)  Nullable(4)
 Nullable(5)  Nullable(6)
 3x2 DataCubes.AbstractArrayWrapper{Nullable{ASCIIString},2,DataCubes.BroadcastAxis{Nullable{ASCIIString},2,DataCubes.AbstractArrayWrapper{Nullable{ASCIIString},1,Array{Nullable{ASCIIString},1}},DataCubes.DictArray{Symbol,2,DataCubes.AbstractArrayWrapper{T,2,A<:AbstractArray{T,N}},Nullable{T}}}}:
 Nullable("X")  Nullable("X")
 Nullable("Y")  Nullable("Y")
 Nullable("Z")  Nullable("Z")

julia> pick(t, [:a, :k])
3 x 2 DictArray

a k |a k
----+----
1 X |2 X
3 Y |4 Y
5 Z |6 Z
```

"""
function pick end

pick(arr::DictArray, field) = selectfield(arr, field)
pick(arr::DictArray, fields::Tuple) = [selectfield(arr, field) for field in fields]
pick(arr::DictArray, field1, fields...) = [selectfield(arr, field) for field in [field1;fields...]]
pick(arr::DictArray, fields::AbstractArray) = selectfields(arr, fields...)
pick(arr::LabeledArray, field) = selectfield(arr, field)
pick(arr::LabeledArray, fields::Tuple) = [selectfield(arr, field) for field in fields]
pick(arr::LabeledArray, fields...) = [selectfield(arr, field) for field in fields]
pick(arr::LabeledArray, fields::AbstractArray) = selectfields(arr, fields...)
pick(arr::LDict, ks::AbstractArray) = selectkeys(arr, ks...)
pick(arr::LDict, k) = arr[k]
pick(arr::LDict, ks::Tuple) = [arr[k] for k in ks]
pick(arr::LDict, ks...) = wraparray_if_nullable_abstractarray(pick(arr, ks))
wraparray_if_nullable_abstractarray{T<:Nullable}(arr::AbstractArray{T}) = wrap_array(arr)
wraparray_if_nullable_abstractarray(arr) = arr

"""

Pick axes from a `LabeledArray`.

* `pickaxis(arr::LabeledArray)` picks a tuple of axes of `arr`.
* `pickaxis(arr::LabeledArray, index::Integer)` picks the axis along the `index` direction.
* `pickaxis(arr::LabeledArray, index::Integer, fields...)` picks the fields `fields` from the `index`th axis in `arr`. It is the same as `pick(pickaxis(arr, index), fields...)`.

##### Examples

```julia
julia> t = larr(a=[1 2;3 4;5 6], b=[:x :y;:z :u;:v :w], axis1=darr(k=["X","Y","Z"]), axis2=[:A,:B])
3 x 2 LabeledArray

  |A   |B
--+----+----
k |a b |a b
--+----+----
X |1 x |2 y
Y |3 z |4 u
Z |5 v |6 w

julia> pickaxis(t)
(3 DictArray

k
--
X
Y
Z
,[Nullable(:A),Nullable(:B)])

julia> pickaxis(t, 1)
3 DictArray

k
--
X
Y
Z

julia> pickaxis(t, 1, :k)
3-element DataCubes.AbstractArrayWrapper{Nullable{ASCIIString},1,Array{Nullable{ASCIIString},1}}:
 Nullable("X")
 Nullable("Y")
 Nullable("Z")
```
"""
function pickaxis end

pickaxis(arr::LabeledArray) = arr.axes
pickaxis(arr::LabeledArray, index::Integer) = arr.axes[index]
pickaxis(arr::LabeledArray, index::Integer, fields...) = pick(arr.axes[index], fields...)

"""

Delete keys or fields from `LDict`, `DictArray`, or `LabeledArray`.

* `delete(dict::LDict, keys...)` deletes the keys `keys` from `dict` and returns a new `LDict`.
* `delete(arr::DictArray, fieldnames...)` deletes the fields for `fieldnames` from `arr` and returns a new `DictArray`.
* `delete(arr::LabeledArray, fieldnames...)` deletes the fields for `fieldnames` from `arr`, either from the base or axes, and returns a new `LabeledArray`.

##### Examples

```julia
julia> delete(LDict(:a=>1, :b=>2, :c=>3), :a, :c)
DataCubes.LDict{Symbol,Int64} with 1 entry:
  :b => 2

julia> delete(darr(a=[1,2,3], b=[:m,:n,:p]), :b)
3 DictArray

a
--
1
2
3

julia> t = larr(a=[1 2;3 4;5 6], b=[:x :y;:z :u;:v :w], axis1=darr(k=["X","Y","Z"]), axis2=[:A,:B])
3 x 2 LabeledArray

  |A   |B
--+----+----
k |a b |a b
--+----+----
X |1 x |2 y
Y |3 z |4 u
Z |5 v |6 w

julia> delete(t, :k, :b)
3 x 2 LabeledArray

  |A |B
--+--+--
  |a |a
--+--+--
1 |1 |2
2 |3 |4
3 |5 |6
```

"""
function delete end

delete(dict::LDict, keys...) = deletekeys(dict, keys...)
delete(arr::DictArray, fieldnames...) = deletefields(arr, fieldnames...)
delete(arr::LabeledArray, fieldnames...) = begin
  newdata = delete(peel(arr), fieldnames...)
  newaxes = map(axis -> delete_fieldnames_from_axis(axis, fieldnames), pickaxis(arr))
  LabeledArray(newdata, newaxes)
end
delete(arr::AbstractArray, fieldnames...) = arr

# a helper function for delete.
delete_fieldnames_from_axis(axis::DictArray, ks) = begin
  newdict = deletekeys(peel(axis), ks...)
  if isempty(newdict)
    DefaultAxis(length(axis))
  else
    DictArray(newdict)
  end
end
delete_fieldnames_from_axis(axis::AbstractArray, ks) = axis

@generated dropna{T,N}(arr::AbstractArrayWrapper{Nullable{T},N}) = quote
  coords_toshow = [falses(n) for n in size(arr)]
  @nloops $N i arr begin
    isna = @nref($N,arr,i).isnull
    @nexprs $N j->(coords_toshow[j][i_j] |= !isna)
  end
  arr[coords_toshow...]
end

@generated dropna_coords{N}(arr::DictArray{TypeVar(:T),N}) = quote
  coords_toshow = [falses(n) for n in size(arr)]
  @nloops $N i arr begin
    isna = all(map(x->x.isnull, getindexvalue(arr, @ntuple($N,i)...)))
    @nexprs $N j->(coords_toshow[j][i_j] |= !isna)
  end
  coords_toshow
end

"""

Remove any `NA` entries. If all elements are `NA` along some slice, that slice will be removed and the array size will shrink.

##### Examples

```julia
julia> t = @darr(a=[1 2 NA;NA 5 NA], b=[NA :n NA;:x NA NA])
2 x 3 DictArray

a b |a b |a b
----+----+----
1   |2 n |
  x |5   |

julia> dropna(t)
2 x 2 DictArray

a b |a b
----+----
1   |2 n
  x |5

julia> m = @larr(a=[1 2 NA;NA 5 NA], b=[NA :n NA;:x NA NA], axis1[:M,:N])
d2 x 3 LabeledArray

  |1   |2   |3
--+----+----+----
  |a b |a b |a b
--+----+----+----
M |1   |2 n |
N |  x |5   |

julia> dropna(m)
2 x 2 LabeledArray

  |1   |2
--+----+----
  |a b |a b
--+----+----
M |1   |2 n
N |  x |5
```

"""
function dropna end

dropna{N}(arr::DictArray{TypeVar(:T),N}) = arr[dropna_coords(arr)...]
dropna{N}(arr::LabeledArray{TypeVar(:T),N}) = arr[dropna_coords(arr.data)...]

"""

`mapna(f::Function, args...)`

Apply `f` to the nullable arrays `args`. It works similarly as `map(f, args...)` but unwraps Nullable from `args`. If any of elements are Nullable, `f` is Nullable, too.

##### Arguments

* `f::Function`: a function to apply.
* `args`: nullable arrays.

##### Returns

A nullable array after applying `f` to elements of `args` for each index. `f` maps non-nullable value to either non-nullable or nullable one. If mapped to a non-nullable value, it will be wrapped by `Nullable` implicitly. If any element of `args` is `NA`, then the return value at that position will be `NA`, too.

##### Examples

```julia
julia> mapna((x,y)->x+y+1, @nalift([1 2 3;4 5 NA]), @nalift([NA 2 3;4 NA NA]))
2x3 DataCubes.AbstractArrayWrapper{Nullable{Int64},2,Array{Nullable{Int64},2}}:
 Nullable{Int64}()  Nullable(5)        Nullable(7)
 Nullable(9)        Nullable{Int64}()  Nullable{Int64}()

julia> mapna((x,y)->Nullable(x+y+1), @nalift([1 2 3;4 5 NA]), @nalift([NA 2 3;4 NA NA]))
2x3 DataCubes.AbstractArrayWrapper{Nullable{Int64},2,Array{Nullable{Int64},2}}:
 Nullable{Int64}()  Nullable(5)        Nullable(7)
 Nullable(9)        Nullable{Int64}()  Nullable{Int64}()
```

"""
function mapna end

mapna(f::Function, args...) = mapna_inner(f, mapna_typecheck(f, args...), args...)

mapna_typecheck(f::Function, args...) = begin
  for xs in zip(args...)
    if all(x->!x.isnull, xs)
      #typ = typeof(f(map(x->x.value, xs)...))
      #return typ <: Nullable ? eltype(typ) : typ
      return mapna_get_eltype_nullable(f(map(x->x.value, xs)...))
    end
  end
  return Any
end

mapna_get_eltype_nullable{T}(::Nullable{T}) = T
mapna_get_eltype_nullable(::Nullable) = Any
mapna_get_eltype_nullable{T}(::T) = T

mapna_inner{T}(f::Function, typ::Type{T}, args...) = map(args...) do xs...
  if any(x->x.isnull, xs)
    Nullable{T}()
  else
    apply_nullable(f(map(x->x.value, xs)...))
  end
end
mapna_typecheck(f::Function, arg) = begin
  for x in arg
    if !x.isnull
      return typeof(f(x.value))
    end
  end
  return Any
end
mapna_inner{T}(f::Function, typ::Type{T}, arg) = map(arg) do x
  if x.isnull
    Nullable{T}()
  else
    apply_nullable(f(x.value))
  end
end
apply_nullable(x::Nullable) = x
apply_nullable(x::DictArray) = x
apply_nullable(x::LabeledArray) = x
apply_nullable(x) = Nullable(x)

"""

`tensorprod(arrs...)`

Calculate the tensor product of `arrs`.

* `tensorprod(::AbstractArray...)` calculates the tensor product where the product operation creates a tuple of the input elements at each position, if every input element at that position is not `NA`. Otherwise, the result will be `NA`.
* `tensorprod(::DictArray...)` calculates the tensor product where the product operation creates a merged `LDict` of the input `LDict`s at each position.
* `tensorprod(::LabeledArray...)` calculates the tensor product of the bases of the inputs, which will be used as the base of the return `LabeledArray`. The axes of the return value will be appropriately extended version of the inputs.

##### Examples

```julia
julia> tensorprod(@nalift([1,2,NA]), @nalift([3,NA]))
3x2 DataCubes.AbstractArrayWrapper{Nullable{Tuple{Int64,Int64}},2,Array{Nullable{Tuple{Int64,Int64}},2}}:
 Nullable((1,3))                 Nullable{Tuple{Int64,Int64}}()
 Nullable((2,3))                 Nullable{Tuple{Int64,Int64}}()
 Nullable{Tuple{Int64,Int64}}()  Nullable{Tuple{Int64,Int64}}()

julia> tensorprod(@darr(a=[1,2,NA]), @darr(b=[3,NA]))
3 x 2 DictArray

a b |a b
----+----
1 3 |1
2 3 |2
  3 |

julia> tensorprod(@larr(a=[1,2,NA], axis1[:m,:n,:p]), @larr(b=[3,NA], axis1[:X,:Y]))
3 x 2 LabeledArray

  |X   |Y
--+----+----
  |a b |a b
--+----+----
m |1 3 |1
n |2 3 |2
p |  3 |
```

"""
function tensorprod end

tensorprod() = error("provide at least one AbstractArray argument.")
tensorprod(arrs::DictArray...) = arrayadd(tensorprod_inner(arrs...)...)
tensorprod(arrs::LabeledArray...) = begin
  newdata = tensorprod(map(x->x.data, arrs)...)
  newaxes = [map(x->collect(x.axes), arrs)...;]
  LabeledArray(newdata, (newaxes...))
end
tensorprod{T<:Nullable}(arrs::AbstractArray{T}...) = begin
  lifted_arrs = tensorprod_inner(arrs...)
  restype = Tuple{map(x->eltype(eltype(x)), lifted_arrs)...}
  map(lifted_arrs...) do elems...
    if any(x->x.isnull, elems)
      Nullable{restype}()
    else
      Nullable(map(x->x.value, elems))
    end
  end
end
tensorprod_inner(arrs::AbstractArray...) = begin
  full_dims = [map(x->collect(size(x)), arrs)...;]
  ndims_vec = [ndims(x) for x in arrs]
  cum_ndims_vec = circshift(cumsum(ndims_vec), 1)
  cum_ndims_vec[1] = 0
  lifted_arrs = map(1:length(arrs)) do i
    arr1d = reshape(arrs[i], length(arrs[i]))
    ndimsi = ndims_vec[i]
    inner_dims = prod(full_dims[1:cum_ndims_vec[i]])
    outer_dims = prod(full_dims[1+cum_ndims_vec[i]+ndimsi:end])
    reshape(repeat(arr1d, inner=[inner_dims], outer=[outer_dims]), full_dims...)
  end
  lifted_arrs
end

# internal use for extract/discard.
immutable RepeatRange
  n::Int
end

immutable Indices{V}
  inds::V
end

"""

`gtake(arr, N1, N2, ...)`

Take a block of an array. similar to `take` in one dimensional case, but is slightly different and more general.
It can also be applied to an `LDict`.
It takes first `N1` elements along direction 1, and similarly for other directions. Repeats if the number of elements are less than `N*`. Picks from rear if `N*` is negative.

##### Examples

```julia
julia> t = larr(a=rand(5,3), b=reshape(1:15,5,3), axis1=[:X,:Y,:Z,:U,:V])
5 x 3 LabeledArray

  |1                     |2                      |3
--+----------------------+-----------------------+-----------------------
  |a                   b |a                   b  |a                   b
--+----------------------+-----------------------+-----------------------
X |0.3219487839233375  1 |0.4863723989946185  6  |0.8784616074632225  11
Y |0.04069063166302023 2 |0.06614308437642014 7  |0.31870618693881947 12
Z |0.7855545407740521  3 |0.5208010912357377  8  |0.4421485355996708  13
U |0.8134241459627629  4 |0.8256022894268482  9  |0.3127049127123851  14
V |0.8536688845922342  5 |0.7263660648355621  10 |0.9315379228053462  15

julia> gtake(t, 3, 2)
3 x 2 LabeledArray

  |1                     |2
--+----------------------+----------------------
  |a                   b |a                   b
--+----------------------+----------------------
X |0.3219487839233375  1 |0.4863723989946185  6
Y |0.04069063166302023 2 |0.06614308437642014 7
Z |0.7855545407740521  3 |0.5208010912357377  8

julia> gtake(t, 3, 4)
3 x 4 LabeledArray

  |1                     |2                     |3                      |4
--+----------------------+----------------------+-----------------------+----------------------
  |a                   b |a                   b |a                   b  |a                   b
--+----------------------+----------------------+-----------------------+----------------------
X |0.3219487839233375  1 |0.4863723989946185  6 |0.8784616074632225  11 |0.3219487839233375  1
Y |0.04069063166302023 2 |0.06614308437642014 7 |0.31870618693881947 12 |0.04069063166302023 2
Z |0.7855545407740521  3 |0.5208010912357377  8 |0.4421485355996708  13 |0.7855545407740521  3

julia> gtake(t, -2, -1)
2 x 1 LabeledArray

  |1
--+----------------------
  |a                  b
--+----------------------
U |0.3127049127123851 14
V |0.9315379228053462 15
```

"""
function gtake end

gtake(arr::LabeledArray, ns...) = extract_inner(arr, map(n->isa(n,Integer) ? RepeatRange(n) : n, ns))
gtake(arr::AbstractArray, ns...) = peel(gtake(larr(arr), ns...))
gtake{K,V}(arr::LDict{K,Nullable{V}}, n::Integer) = gtake_integer(arr, n)
gtake(arr::LDict, n::Integer) = gtake_integer(arr, n)
gtake_integer(arr::LDict, n::Integer) = begin
  lenarr = length(arr)
  if abs(n) > lenarr
    throw(ArgumentError("LDict length $(lenarr) is smaller than the range $(n)"))
  end
  if n >= 0
    LDict(arr.keys[1:n], arr.values[1:n])
  else
    LDict(arr.keys[lenarr+n+1:end], arr.values[lenarr+n+1:end])
  end
end

"""

`extract(arr, ns...)`

Extract a block of array using labels for `LabeledArray`s, indices for other types of arrays or `LDict`s.

##### Arguments

* `arr`: an `AbstractArray` or `LDict`.
* `ns...`: each element in `ns` chooses specific elements along that direction. The element can be
  * `Colon()` (`:`): the entire range.
  * a label along that direction. If the axis along the direction is `DictArray`, the label can be either an `LDict` for its element or a tuple to denote the values of `LDict`.
  * array of labels along that direction.
  * a boolean array of the same size as the axis along that direction to denote which position to choose.
  * a function that takes the axis along that direction and generates either an array of integers or a boolean array for the selected positions.

##### Return

An array or `LDict` of the same type as `arr`, which is selected based on `ns`.... All indices will be chosen for the rest of the directions not specified in `ns`.... If any label is missing or the integer range is out of bound, `NA` will be used for that element in the return value. If an element in `ns` is scalar, the dimension along that direction will be collapsed just as in `slice`.

##### Examples

```julia
julia> t = larr(a=map(x->'a'+x,reshape(0:14,5,3)), b=reshape(1:15,5,3), axis1=[:X,:Y,:Z,:U,:V], axis2=darr(r1=[:A,:A,:B],r2=[:m,:n,:n]))
5 x 3 LabeledArray

r1 |A   |A    |B
r2 |m   |n    |n
---+----+-----+-----
   |a b |a b  |a b
---+----+-----+-----
X  |a 1 |f 6  |k 11
Y  |b 2 |g 7  |l 12
Z  |c 3 |h 8  |m 13
U  |d 4 |i 9  |n 14
V  |e 5 |j 10 |o 15

julia> extract(t, [:X,:V,:W], map(Nullable,(:A,:m)))
3 LabeledArray

  |a b
--+----
X |a 1
V |e 5
W |

julia> extract(t, [:X,:V,:W], darr(r1=[:A,:B],r2=[:m,:m]))
3 x 2 LabeledArray

r1 |A   |B
r2 |m   |m
---+----+----
   |a b |a b
---+----+----
X  |a 1 |
V  |e 5 |
W  |    |

julia> extract(t, :, d->d[:r1] .== :A)
5 x 2 LabeledArray

r1 |A   |A
r2 |m   |n
---+----+-----
   |a b |a b
---+----+-----
X  |a 1 |f 6
Y  |b 2 |g 7
Z  |c 3 |h 8
U  |d 4 |i 9
V  |e 5 |j 10
```

"""
function extract end

extract(arr::LabeledArray, ns...) = extract_inner(arr, ns)
extract(arr::AbstractArray, ns...) = peel(extract(larr(arr), ns...))
extract{K,V}(ldict::LDict{K,Nullable{V}}, k) = begin
  ind = findfirst(ldict.keys, k)
  ind > 0 ? ldict.values[ind] : Nullable{V}()
end
extract{K}(ldict::LDict{K}, k) = begin
  ind = findfirst(ldict.keys, k)
  ind > 0 ? ldict.values[ind] : Nullable{Any}()
end
extract{K,V}(ldict::LDict{K,Nullable{V}}, ks::AbstractArray) = begin
  ldictkeys = ldict.keys
  ldictvalues = ldict.values
  LDict(collect(ks), [(ind=findfirst(ldictkeys,k);ind>0 ? ldictvalues[ind] : Nullable{V}()) for k in ks])
end
extract{K}(ldict::LDict{K}, ks::AbstractArray) = begin
  ldictkeys = ldict.keys
  ldictvalues = ldict.values
  LDict(collect(ks), [(ind=findfirst(ldictkeys,k);ind>0 ? ldictvalues[ind] : Nullable{Any}()) for k in ks])
end
extract(ldict::LDict{Function,Nullable}, ks::Function) = error("not yet implemented.")
# to disambiguate a method choice.
extract{K,V}(ldict::LDict{K,Nullable{V}}, ks::Function) = extract_function_inner(ldict, ks)
extract(ldict::LDict, ks::Function) = extract_function_inner(ldict, ks)
extract_function_inner(ldict::LDict, ks::Function) = begin
  # it does not make sense to automatically extend the keys, because they are not Nullable.
  inds = convert_boolarray_if_necessary(ks(keys(ldict)))
  extract_ldict_from_indices(ldict, inds)
end
extract_ldict_from_indices(ldict::LDict, inds::AbstractArray) = LDict(keys(ldict)[inds], values(ldict)[inds])
extract_ldict_from_indices(ldict::LDict, ind) = values(ldict)[ind]

extract_inner(arr::LabeledArray, ns::Tuple) = begin
  lenns = length(ns)
  extended_ns_tuples = ntuple(ndims(arr)) do d
    if d <= lenns
      process_one_range(arr.axes[d], ns[d])
    else
      (Colon(), Colon())
    end
  end
  #scalar_directions_type = ([(isa(v, AbstractArray) || isa(v, Function) || isa(v, Colon) || isa(v, RepeatRange)) ? Colon() : 0 for v in extended_ns]...)
  extended_ns = map(x->x[1], extended_ns_tuples)
  scalar_directions_type = map(x->x[2], extended_ns_tuples)
  indices = [create_extract_indices(axis, region) for (axis,region) in zip(arr.axes, extended_ns)]
  choose_array(arr, (indices...), scalar_directions_type)
end

process_one_range(axis, x::AbstractArray) = (x, Colon())
process_one_range(axis, x::Colon) = (x, Colon())
process_one_range(axis, x::RepeatRange) = (x, Colon())
process_one_range(axis, x) = (x, 0)
process_one_range(axis, x::Function) = process_one_range_indices(convert_boolarray_if_necessary(x(axis)))
process_one_range_indices(x::AbstractArray) = (Indices(x), Colon())
process_one_range_indices(x) = (Indices([x]), 0)

"""

`gdrop(arr, N1, N2, ...)`

Drop a block of an array. similar to `drop` in one dimensional case, but is slightly different and more general.
It can also be applied to an `LDict`.
It drops the first `N1` elements along direction 1, and similarly for other directions. Drops from rear if `N*` is negative.

##### Examples

```julia
julia> t = larr(a=rand(5,3), b=reshape(1:15,5,3), axis1=[:X,:Y,:Z,:U,:V])
5 x 3 LabeledArray

  |1                     |2                      |3
--+----------------------+-----------------------+-----------------------
  |a                   b |a                   b  |a                   b
--+----------------------+-----------------------+-----------------------
X |0.27289790581491746 1 |0.8493197848353495  6  |0.8370920536703472  11
Y |0.8424940964507834  2 |0.21518951524950136 7  |0.9290437789813346  12
Z |0.9498541774517255  3 |0.942687447396005   8  |0.1341678643795654  13
U |0.7356663426240728  4 |0.7662948222160162  9  |0.24109069576951692 14
V |0.8716491751450759  5 |0.27472373001295436 10 |0.08909928028262804 15

julia> gdrop(t, 3, 2)
2 x 1 LabeledArray

  |1
--+-----------------------
  |a                   b
--+-----------------------
U |0.24109069576951692 14
V |0.08909928028262804 15

julia> gdrop(t, 5)
0 x 3 LabeledArray

 |1   |2   |3
-+----+----+----
 |a b |a b |a b

julia> gdrop(t, -3, -2)
2 x 1 LabeledArray

  |1
--+----------------------
  |a                   b
--+----------------------
X |0.27289790581491746 1
Y |0.8424940964507834  2
```

"""
function gdrop end

gdrop(arr::LabeledArray, ns...) = discard_inner(arr, map(n->isa(n,Integer) ? RepeatRange(n) : n, ns))
gdrop(arr::AbstractArray, ns...) = peel(gdrop(larr(arr), ns...))

gdrop{K,V}(arr::LDict{K,Nullable{V}}, n::Integer) = gdrop_integer(arr, n)
gdrop(arr::LDict, n::Integer) = gdrop_integer(arr, n)
gdrop_integer(arr::LDict, n::Integer) = begin
  if n >= 0
    LDict(arr.keys[n+1:end], arr.values[n+1:end])
  else
    LDict(arr.keys[1:end+n], arr.values[1:end+n])
  end
end

"""

`discard(arr, ns...)`

Discard a block of array discarding all elements specified by `ns`..., using labels for `LabeledArray`s, indices for other types of arrays or `LDict`s.

##### Arguments

* `arr`: an `AbstractArray` or `LDict`.
* `ns...`: each element in `ns` chooses specific elements along that direction. The element can be
  * `Colon()` (`:`): the entire range will be removed and the return value will be empty.
  * a label or array of labels along that direction to discard.
  * a boolean array of the same size as the axis along that direction to denote which position to discard.
  * a function that takes the axis along that direction and generates either an array of integers or a boolean array for the deleted positions.

##### Return

An array or `LDict` of the same type as `arr`, which is selected based on `ns`.... All indices will be chosen for the rest of the directions not specified in `ns`.... If any label is missing or the integer range is out of bound, it will be ignored.

##### Examples

```julia
julia> t = larr(a=map(x->'a'+x,reshape(0:14,5,3)), b=reshape(1:15,5,3), axis1=[:X,:Y,:Z,:U,:V], axis2=darr(r1=[:A,:A,:B],r2=[:m,:n,:n]))
5 x 3 LabeledArray

r1 |A   |A    |B
r2 |m   |n    |n
---+----+-----+-----
   |a b |a b  |a b
---+----+-----+-----
X  |a 1 |f 6  |k 11
Y  |b 2 |g 7  |l 12
Z  |c 3 |h 8  |m 13
U  |d 4 |i 9  |n 14
V  |e 5 |j 10 |o 15

julia> discard(t, [:X,:V,:W], map(Nullable,(:A,:m)))
3 x 2 LabeledArray

r1 |A   |B
r2 |n   |n
---+----+-----
   |a b |a b
---+----+-----
Y  |g 7 |l 12
Z  |h 8 |m 13
U  |i 9 |n 14

julia> discard(t, [:X,:V,:W], darr(r1=[:A,:B],r2=[:m,:m]))
3 x 2 LabeledArray

r1 |A   |B
r2 |n   |n
---+----+-----
   |a b |a b
---+----+-----
Y  |g 7 |l 12
Z  |h 8 |m 13
U  |i 9 |n 14

julia> discard(t, [], d->d[:r1] .== :A)
5 x 1 LabeledArray

r1 |B
r2 |n
---+-----
   |a b
---+-----
X  |k 11
Y  |l 12
Z  |m 13
U  |n 14
V  |o 15
```

"""
function discard end

discard(arr::LabeledArray, ns...) = discard_inner(arr, ns)
discard(arr::AbstractArray, ns...) = peel(discard(larr(arr), ns...))

discard(ldict::LDict, k) = deletekeys(ldict, k)
discard(ldict::LDict, ks::AbstractArray) = deletekeys(ldict, ks...)
discard(ldict::LDict, ks::Function) = begin
  inds = convert_boolarray_if_necessary(ks(keys(ldict)))
  lendict = length(ldict)
  indices = trues(lendict)
  for i in inds
    if checkbounds(Bool, lendict, i)
      indices[i] = false
    end
  end
  LDict(keys(ldict)[indices], values(ldict)[indices])
end

discard_inner(arr::LabeledArray, ns::Tuple) = begin
  lenns = length(ns)
  extended_ns = ntuple(ndims(arr)) do d
    if d <= lenns
      ns[d]
    else
      []
    end
  end
  indices = [create_discard_indices(axis, region) for (axis,region) in zip(arr.axes, extended_ns)]
  choose_array(arr, (indices...), ntuple(d->Colon(), ndims(arr)))
end

choose_array{N}(arr::LabeledArray{TypeVar(:T),N}, indices::NTuple{N}, scalar_directions_type) = begin
  non_scalar_directions = find([x==Colon() for x in scalar_directions_type])
  locs_tuple = map(x->x[2], indices)
  newdata = choose_array(peel(arr), locs_tuple, scalar_directions_type)
  newaxis = ([x[1] for x in indices][non_scalar_directions]...)
  if isempty(non_scalar_directions)
    newdata
  else
    LabeledArray(newdata, newaxis)
  end
end

choose_array{N}(arr::DictArray{TypeVar(:K),N}, indices::NTuple{N}, scalar_directions_type) = begin
  non_scalar_directions = find([x==Colon() for x in scalar_directions_type])
  if isempty(non_scalar_directions)
    mapvalues(v->choose_array(v, indices, scalar_directions_type), peel(arr))
  else
    DictArray(mapvalues(v->choose_array(v, indices, scalar_directions_type), peel(arr)))
  end
end

@generated choose_array{T,N}(arr::AbstractArray{Nullable{T},N}, indices::NTuple{N}, scalar_directions_type) = begin
  scalar_directions = find([t == Int for t in scalar_directions_type.types])
  non_scalar_directions = find([t == Colon for t in scalar_directions_type.types])
  if isempty(non_scalar_directions)
    # the result is a scalar, not an array.
    quote
      coords = @ntuple($N, d->indices[d][1])
      any(x->x==0, coords) ? Nullable{T}() : arr[coords...]
    end
  else
    result_coords_tuple = Expr(:tuple, [symbol("i_", d) for d in non_scalar_directions]...)
    quote
      resultsize = [length(ind) for ind in indices][$non_scalar_directions]
      result = similar(arr, resultsize...)
      @nloops $N i d->eachindex(indices[d]) begin
        coords = @ntuple($N,d->indices[d][i_d])
        result[$result_coords_tuple...] = any(x->x==0, coords) ? Nullable{T}() : arr[coords...]
      end
      result
    end
  end
end

# if the region is a positive integer, take the first region elements.
# if there are less than region elements, repeat the whole array repeatedly.
# if negative, take the last region elements. repeat if necessary.
create_extract_indices(axis::AbstractVector, region0::RepeatRange) = begin
  region = region0.n
  lenaxis = length(axis)
  indices::Array{Int,1}
  result = if region >= 0
    if region <= lenaxis
      1:region
    else
      indices = Array(Int, region)
      pos = 1
      for i in eachindex(indices)
        indices[i] = pos
        pos += 1
        if pos > lenaxis
          pos = 1
        end
      end
      indices
    end
  else
    nregion = -region
    if nregion <= lenaxis
      lenaxis-nregion+1:lenaxis
    else
      indices = Array(Int, nregion)
      pos = lenaxis
      for i in length(indices):-1:1
        indices[i] = pos
        pos -= 1
        if pos == 0
          pos = lenaxis
        end
      end
      indices
    end
  end
  (isa(axis, DefaultAxis) ? DefaultAxis(length(result)) : axis[result]), result
end
create_extract_indices(axis::AbstractVector{Nullable{RepeatRange}}, region0::RepeatRange) = error("this should not be called")
create_extract_indices(axis::AbstractVector{Nullable{Colon}}, region::Colon) = axis, eachindex(axis)
create_extract_indices(axis::AbstractVector, region::Colon) = axis, eachindex(axis)
create_extract_indices(axis::DefaultAxis, region::AbstractArray{Nullable{Int}}) = begin
  lenaxis = length(axis)
  indices = map(region) do i
    if i.isnull
      0
    elseif checkbounds(Bool, lenaxis, i.value)
      i.value
    else
      0
    end
  end
  DefaultAxis(length(indices)), indices
end
create_extract_indices(axis::DefaultAxis, region::AbstractArray{Int}) = begin
  lenaxis = length(axis)
  if checkbounds(Bool, lenaxis, region)
    DefaultAxis(length(region)), region
  else
    indices = map(region) do i
      if checkbounds(Bool, lenaxis, i)
        i
      else
        0
      end
    end
    DefaultAxis(length(indices)), indices
  end
end
create_extract_indices(axis::DictArray, region::LDict) = create_extract_indices(axis, DictArray(mapvalues(v->wrap_array([v]), region)))
create_extract_indices(axis::DictArray, region::Tuple) = create_extract_indices(axis, DictArray(keys(axis), [wrap_array([v]) for v in region]))
create_extract_indices(axis::DictArray, region::DictArray) =
  create_extract_indices_typed(Tuple{[elvaluetype(v) for v in values(axis)]...}, axis, region)
create_extract_indices_typed{T}(::Type{T}, axis::DictArray, region::DictArray) = begin
  @assert(keys(axis) == keys(region))
  lenaxis = length(axis)
  indmap = Dict{T,Int}()
  for i in 1:lenaxis
    value::T = getindexvalue(axis, T, i)
    indmap[value] = i
  end
  indices = map(1:length(region)) do i
    value::T = getindexvalue(region, T, i)
    if haskey(indmap, value)
      indmap[value]
    else
      0
    end
  end
  region, indices
end
create_extract_indices{T}(axis::AbstractVector{Nullable{T}}, region::AbstractArray{T}) = create_extract_indices(axis, nalift(region))
create_extract_indices{T}(axis::AbstractVector{Nullable{T}}, region::Nullable{T}) = create_extract_indices(axis, wrap_array([region]))
create_extract_indices{T}(axis::AbstractVector{Nullable{T}}, region::T) = create_extract_indices(axis, wrap_array([Nullable(region)]))
create_extract_indices{T<:Indices}(axis::AbstractVector{Nullable{T}}, region::T) = error("this should not happen.")
#create_extract_indices(axis::DefaultAxis, region::Function) = create_extract_indices(axis, convert_boolarray_if_necessary(region(axis)))
create_extract_indices(axis::DefaultAxis, region::Indices) = create_extract_indices(axis, region.inds)
create_extract_indices(axis::DictArray, region::Indices) = create_extract_indices_typed(axis, region, Tuple{[elvaluetype(v) for v in values(axis)]...})
create_extract_indices_typed{T}(axis::DictArray, region::Indices, ::Type{T}) = begin
  converted_indices = collect(region.inds)
  converted_region = similar(axis, length(converted_indices))
  lenaxis = length(axis)
  for i in eachindex(converted_indices)
    if checkbounds(Bool, lenaxis, converted_indices[i])
      converted_region[i] = getindexvalue(axis, T, converted_indices[i])
    else
      setna!(converted_region, i)
      converted_indices[i] = 0
    end
  end
  converted_region, converted_indices
end

create_extract_indices{T}(axis::AbstractVector{Nullable{T}}, region::Indices) = begin
  converted_indices = collect(region.inds)
  converted_region = similar(axis, length(converted_indices))
  lenaxis = length(axis)
  for i in eachindex(converted_indices)
    if checkbounds(Bool, lenaxis, converted_indices[i])
      converted_region[i] = axis[converted_indices[i]]
    else
      converted_region[i] = Nullable{T}()
      converted_indices[i] = 0
    end
  end
  converted_region, converted_indices
end
convert_boolarray_if_necessary(indices::AbstractArray{Bool,1}) = find(indices)
convert_boolarray_if_necessary(indices::AbstractArray{Nullable{Bool},1}) = find(ignabool(indices))
convert_boolarray_if_necessary(indices) = indices

create_extract_indices{T}(axis::AbstractVector{T}, region::AbstractArray{T}) = begin
  lenaxis = length(axis)
  indmap = Dict{T,Int}()
  for i in 1:lenaxis
    value::T = axis[i]
    indmap[value] = i
  end
  indices = map(1:length(region)) do i
    value::T = region[i]
    if haskey(indmap, value)
      indmap[value]
    else
      0
    end
  end
  region, indices
end

# if the region is a positive integer, drop the first region elements.
# if negative, drop the last region elements.
create_discard_indices(axis::AbstractArray, region0::RepeatRange) = begin
  region = region0.n
  lenaxis = length(axis)
  indices = if region >= 0
    region+1:lenaxis
  else
    1:lenaxis+region
  end
  axis[indices], indices
end
create_discard_indices(axis::DefaultAxis, region0::RepeatRange) = begin
  region = region0.n
  lenaxis = length(axis)
  indices = if region >= 0
    region+1:lenaxis
  else
    1:lenaxis+region
  end
  DefaultAxis(length(indices)), indices
end
create_discard_indices(axis::AbstractVector{Nullable{RepeatRange}}, region0::RepeatRange) = error("this should not be called")
create_discard_indices(axis::AbstractVector{Nullable{Colon}}, region::Colon) = axis[Int[]], Int[]
create_discard_indices(axis::AbstractArray, region::Colon) = axis[Int[]], Int[]
create_discard_indices(axis::DefaultAxis, region::AbstractArray{Nullable{Int}}) = create_discard_indices(axis, Int[x.isnull ? 0 : x.value for x in region])
create_discard_indices(axis::DefaultAxis, region::AbstractArray{Int}) = begin
  if isempty(region)
    return axis, eachindex(axis)
  end
  lenaxis = length(axis)
  indices = trues(lenaxis)
  for i in region
    if checkbounds(Bool, lenaxis, i)
      indices[i] = false
    end
  end
  locs = find(indices)
  DefaultAxis(length(locs)), locs
end
create_discard_indices(axis::DictArray, region::LDict) = create_discard_indices(axis, DictArray(mapvalues(v->wrap_array([v]), region)))
create_discard_indices(axis::DictArray, region::Tuple) = create_discard_indices(axis, DictArray(keys(axis), [wrap_array([v]) for v in region]))
create_discard_indices(axis::DictArray, region::DictArray) =
  create_discard_indices_typed(Tuple{[eltype(v) for v in axis.data]...}, axis, region)
create_discard_indices_typed{T}(::Type{T}, axis::DictArray, region::DictArray) = begin
  if isempty(region)
    return axis, eachindex(axis)
  end
  @assert(keys(axis) == keys(region))
  lenaxis = length(axis)
  indmap = Dict{T,Int}()
  for i in 1:lenaxis
    value::T = getindexvalue(axis, T, i)
    indmap[value] = i
  end
  indices = trues(lenaxis)
  for i in 1:length(region)
    value::T = getindexvalue(region, T, i)
    if haskey(indmap, value)
      indices[indmap[value]] = false
    end
  end
  locs = find(indices)
  axis[locs], locs
end
promote_to_array_if_necessary(x::AbstractArray) = x
promote_to_array_if_necessary(x::Colon) = x
promote_to_array_if_necessary(x) = [x]

create_discard_indices(axis::AbstractVector{Nullable{Function}}, region::Function) = error("this is ambiguous")
create_discard_indices{T<:Function}(axis::AbstractVector{Nullable{T}}, region::Function) = error("this is ambiguous")
create_discard_indices{T}(axis::AbstractVector{Nullable{T}}, region::AbstractArray{T}) = create_discard_indices(axis, nalift(region))
create_discard_indices{T}(axis::AbstractVector{Nullable{T}}, region::Nullable{T}) = create_discard_indices(axis, wrap_array([region]))
create_discard_indices{T}(axis::AbstractVector{Nullable{T}}, region::T) = create_discard_indices(axis, wrap_array([Nullable(region)]))
create_discard_indices(axis::DefaultAxis, region::Function) = create_discard_indices(axis, promote_to_array_if_necessary(convert_boolarray_if_necessary(region(axis))))
create_discard_indices(axis::AbstractVector, region::Function) = begin
  converted_indices = collect(promote_to_array_if_necessary(convert_boolarray_if_necessary(region(axis))))
  if isempty(converted_indices)
    return axis, eachindex(axis)
  end
  lenaxis = length(axis)
  indices = trues(lenaxis)
  for i in converted_indices
    if checkbounds(Bool, lenaxis, i)
      indices[i] = false
    end
  end
  locs = find(indices)
  axis[locs], locs
end

create_discard_indices{T}(axis::AbstractVector{T}, region::AbstractArray{T}) = begin
  if isempty(region)
    return axis, eachindex(axis)
  end
  lenaxis = length(axis)
  indmap = Dict{T,Int}()
  for i in 1:lenaxis
    value::T = axis[i]
    indmap[value] = i
  end
  indices = trues(lenaxis)
  for i in 1:length(region)
    value::T = region[i]
    if haskey(indmap, value)
      indices[indmap[value]] = false
    end
  end
  locs = find(indices)
  axis[locs], locs
end
create_discard_indices{T<:AbstractArray}(axis::AbstractVector{Nullable{T}}, region::AbstractArray{T}) = create_discard_indices_generic(axis, region)
create_discard_indices{T<:AbstractArray}(axis::AbstractVector{Nullable{T}}, region::AbstractArray) = create_discard_indices_generic(axis, region)
create_discard_indices{T}(axis::AbstractVector{T}, region::AbstractArray) = create_discard_indices_generic(axis, region)
create_discard_indices_generic{T}(axis::AbstractVector{T}, region::AbstractArray) = if isempty(region)
  create_discard_indices(axis, T[])
else
  throw(ArgumentError("possibly, a wrong type region argument is given: axis $(axis), region $(region)"))
end

"""

`providenames(arr::LabeledArray, create_fieldname::Funcion)`

Add generic field names for fields without field names in a `LabeledArray`.
This makes the data and all the axes components `DictArray`s.
This is useful when you want to apply `selct`/`update`/`leftjoin`/`innerjoin` whose interface is friendlier to `DictArray`s than general `AbstractArray`s.
The reverse operation, removing generic field names, is done by `withdrawnames`.
An optional argument `create_fieldname` is a function that gives a symbol that will be used as a new field name given an integer index.
By default, it generates `:xN` for an index integer `N`.

##### Examples

```julia
julia> t = larr([1 2 3;4 5 6], axis1=[:X,:Y], axis2=darr(k=["A","B","C"]))
2 x 3 LabeledArray

k |A |B |C
--+--+--+--
  |  |  |
--+--+--+--
X |1 |2 |3
Y |4 |5 |6

julia> providenames(t)
2 x 3 LabeledArray

k  |A  |B  |C
---+---+---+---
x2 |x1 |x1 |x1
---+---+---+---
X  |1  |2  |3
Y  |4  |5  |6
```

"""
function providenames end

providenames(arr::LabeledArray, create_fieldname::Function=generic_fieldname) = begin
  tracker = []
  newdata = if isa(arr.data, DictArray)
    arr.data
  else
    newname = create_additional_fieldname(arr, tracker, create_fieldname)
    DictArray(newname => arr.data)
  end
  newaxis = ntuple(length(arr.axes)) do d
    if isa(arr.axes[d], DictArray)
      arr.axes[d]
    else
      newname = create_additional_fieldname(arr, tracker, create_fieldname)
      DictArray(newname => arr.axes[d])
    end
  end
  LabeledArray(newdata, newaxis)
end

"""

`withdrawnames(arr::LabeledArray, check_fieldname::Function)`

Remove generic field names from a `LabeledArray` if possible.

##### Arguments
* `arr` is an input `LabeledArray`.
* `check_fieldname` is a function that returns whether an input field name is a generic name. By default, a field name is generic if it is a symbol of the form `:xN` for some integer `N`.

##### Return
If a field name in `LabeledArray` gives `true` when applied to `check_fieldname`, and the field name can be removed, it is removed in the return `LabeledArray`. A field name can be removed if it is a part of a `DictArray` with only one field.

##### Examples

```julia
julia> t = larr([1 2 3;4 5 6], axis1=[:X,:Y], axis2=darr(k=["A","B","C"]))
2 x 3 LabeledArray

k |A |B |C
--+--+--+--
  |  |  |
--+--+--+--
X |1 |2 |3
Y |4 |5 |6

julia> providenames(t)
2 x 3 LabeledArray

k  |A  |B  |C
---+---+---+---
x2 |x1 |x1 |x1
---+---+---+---
X  |1  |2  |3
Y  |4  |5  |6

julia> withdrawnames(providenames(t))
2 x 3 LabeledArray

k |A |B |C
--+--+--+--
  |  |  |
--+--+--+--
X |1 |2 |3
Y |4 |5 |6
```

"""
function withdrawnames end

withdrawnames(arr::LabeledArray, check_fieldname::Function=is_generic_fieldname) = begin
  newdata = if isa(arr.data, DictArray) && length(arr.data.data) == 1 && check_fieldname(arr.data.data.keys[1])
    arr.data.data.values[1]
  else
    arr.data
  end
  newaxes = ntuple(length(arr.axes)) do d
    arraxesd = arr.axes[d]
    if isa(arraxesd, DictArray) && length(arraxesd.data) == 1 && check_fieldname(arraxesd.data.keys[1])
      arraxesd.data.values[1]
    else
      arraxesd
    end
  end
  LabeledArray(newdata, newaxes)
end

# note the return value of namerge_inner! matters.
namerge_inner!(x::Nullable, y::Nullable) = y.isnull ? x : y
namerge_inner!{T<:Nullable}(x::Nullable, y::AbstractArray{T}) = begin
  result = similar(y)
  for i in eachindex(y)
    result[i] = y[i].isnull ? x : y[i]
  end
  result
end
namerge_inner!{T<:Nullable,U<:Nullable}(x::AbstractArray{T}, y::U) = begin
  if y.isnull
    x
  else
    fill!(x, y)
    x
  end
end
namerge_inner!{T<:Nullable}(x::AbstractArray{T}, y::AbstractArray{T}) = begin
  if size(x) != size(y)
    throw(ArgumentError("array sizes do not match: $(size(x)) and $(size(y))"))
  end
  for i in eachindex(x)
    x[i] = y[i].isnull ? x[i] : y[i]
  end
  x
end
namerge_helper_copy(x::AbstractArray) = copy(x)
namerge_helper_copy(x) = x

namerge() = error("require at least one argument")
namerge(x) = x

"""

`namerge(xs...)`

Combine `Nullable` arrays and `Nullable` elements, having the later arguments override the preceding ones if the new element is not null.

##### Arguments

`xs...` consists of either a `AbstractArrayWrapper` with `Nullable` element type, or a `Nullable` variable. If an element is neither `AbstractArray` or `Nullable`, it will be wrapped by `Nullable`.

##### Return

* When there is no argument, an error will raise.
* If there is only one argument, that argument will be returned.
* If there are two arguments and if the two are not `AbstractArray`, the second argument will be returned only if it is not null. Otherwise, the first argument will be returned.
* If there are two arguments, and the two are `Nullable` arrays, the element at each position will be the one from the first argument if the second argument element is null. Otherwise, the element from the second argument will be used. If any argument is not `AbstractArray`, it will be promoted to a `Nullable` array.

##### Examples

```julia

julia> namerge(10, @nalift([1 2 NA;4 NA NA]))
2x3 DataCubes.AbstractArrayWrapper{Nullable{Int64},2,Array{Nullable{Int64},2}}:
 Nullable(1)  Nullable(2)   Nullable(10)
 Nullable(4)  Nullable(10)  Nullable(10)

julia> namerge(@nalift([1 2 NA;4 NA NA]), 10)
2x3 DataCubes.AbstractArrayWrapper{Nullable{Int64},2,Array{Nullable{Int64},2}}:
 Nullable(10)  Nullable(10)  Nullable(10)
 Nullable(10)  Nullable(10)  Nullable(10)

julia> namerge(10, @nalift([1 2 NA;4 NA NA]))
2x3 DataCubes.AbstractArrayWrapper{Nullable{Int64},2,Array{Nullable{Int64},2}}:
 Nullable(1)  Nullable(2)   Nullable(10)
 Nullable(4)  Nullable(10)  Nullable(10)

julia> namerge(@nalift([1 2 NA;4 NA NA]), @nalift([11 NA NA;14 15 NA]))
2x3 DataCubes.AbstractArrayWrapper{Nullable{Int64},2,Array{Nullable{Int64},2}}:
 Nullable(11)  Nullable(2)   Nullable{Int64}()
 Nullable(14)  Nullable(15)  Nullable{Int64}()

```

"""
function namerge end

namerge(xs...) = begin
  result = nalift(namerge_helper_copy(xs[1]))
  for i in 2:length(xs)
    result = namerge_inner!(result, nalift(xs[i]))
  end
  result
end

elvaluetype(arr::AbstractArray) = elvaluetype(typeof(arr))
elvaluetype{T<:AbstractArray}(::Type{T}) = eltype(T)
elvaluetype{T<:EnumerationArray}(::Type{T}) = Int
elvaluetype{T,N,A}(::Type{FloatNAArray{T,N,A}}) = T
elvaluetype{T,N,A}(::Type{AbstractArrayWrapper{T,N,A}}) = elvaluetype(A)
elvaluetype{T,N,A,I}(::Type{SubArrayView{T,N,A,I}}) = elvaluetype(A)
elvaluetype{T,N,V,W}(::Type{BroadcastAxis{T,N,V,W}}) = elvaluetype(V)

Base.push!{T,V,R}(arr::EnumerationArray{T,1,V,R}, elem::Nullable{T}) = begin
  if elem.isnull
    push!(arr.elems, zero(R))
  else
    # note a T value that is not one of the pool values will be assigned to 0.
    # that is, it is treated as null.
    push!(arr.elems, convert(R, findfirst(arr.pool, elem.value)))
  end
end
Base.push!{T,V,R}(arr::EnumerationArray{T,1,V,R}, elem::R) = push!(arr.elems, elem)
Base.push!{T,A}(arr::FloatNAArray{T,1,A}, elem::Nullable{T}) = push!(arr.data, elem.isnull ? convert(T,NaN) : elem.value)
Base.push!{T,A}(arr::FloatNAArray{T,1,A}, elem::T) = push!(arr.data, elem)


text:
using Git

cd("..")
run(`$gitcmd submodule update --init`)

ext: jl


text:
using Test

using MARY_fRG.Fermi
using MARY_fRG.FlowEquation
using MARY_fRG.Drawers

@testset "测试精细化" begin
    qs = common_triangle_lattice(0.11)
    Γ4 = TFGamma4(qs, 6.0, 24, 10)
    plt = draw_lines(Γ4.model.brillouin.edges)
    draw_polygon!(plt, Γ4.ltris; pcidx=Γ4.lpats)
    @savepng plt "refine1"
    #
    for idx in 2:1:9
        newΓ4 = TFGamma4_refine_ltris_mt(Γ4, 5.0)
        println(idx, ": ", length(newΓ4.ltris))
        plt = draw_lines(newΓ4.model.brillouin.edges)
        draw_polygon!(plt, newΓ4.ltris; pcidx=newΓ4.lpats)
        @savepng plt "refine"*string(idx)
        Γ4 = newΓ4
    end
end


text:
function is_triangle(sides)
  sides = sort(sides)
  sides[1] > 0 && sides[3] < sides[1] + sides[2]
end

function is_equilateral(sides)
  is_triangle(sides) && sides[1] == sides[2] && sides[2] == sides[3]
end

function is_isosceles(sides)
  sides = sort(sides)
  is_triangle(sides) && sides[1] == sides[2] || sides[2] == sides[3]
end

function is_scalene(sides)
  is_triangle(sides) && !is_isosceles(sides)
end

max_stars_repo_name: cmccandless/ExercismSolutions-julia

max_issues_repo_name: cmccandless/ExercismSolutions-julia

max_forks_repo_name: cmccandless/ExercismSolutions-julia


text:
# Note that this script can accept some limited command-line arguments, run
# `julia build_tarballs.jl --help` to see a usage message.
using BinaryBuilder, Pkg

name = "mold"
version = v"1.2.0"

# Collection of sources required to complete build
sources = [
    GitSource("https://github.com/rui314/mold.git", "1f55e40a8a967894816d6366bcc3d08de74c84b4")
]

# Bash recipe for building across all platforms
script = raw"""
cd $WORKSPACE/srcdir/mold/
make -j${nproc}
make install PREFIX="${prefix}" BINDIR="${bindir}" LIBDIR="${libdir}"
"""

# These are the platforms we will build for by default, unless further
# platforms are passed in on the command line
platforms = expand_cxxstring_abis(filter!(Sys.islinux, supported_platforms()))

# The products that we will ensure are always built
products = [
    ExecutableProduct("mold", :mold)
]

# Dependencies that must be installed before this package can be built
dependencies = [
    Dependency(PackageSpec(name="Zlib_jll", uuid="83775a58-1f1d-513f-b197-d71354ab007a"))
    Dependency(PackageSpec(name="OpenSSL_jll", uuid="458c3c95-2e84-50aa-8efc-19380b2a3a95"))
]

# Build the tarballs, and possibly a `build.jl` as well.
build_tarballs(ARGS, name, version, sources, script, platforms, products, dependencies; julia_compat="1.6", preferred_gcc_version=v"10")

max_stars_repo_name: MichelJuillard/Yggdrasil

max_issues_repo_name: MichelJuillard/Yggdrasil

max_forks_repo_name: MichelJuillard/Yggdrasil


text:
from_default_theme(attr) = Makie.current_default_theme()[attr]

function legend_elements(::Type{Scatter};
                         marker=from_default_theme(:marker),
                         markerpoints=[Point2f(0.5, 0.5)],
                         color=from_default_theme(:markercolor),
                         kwargs...)
    return [MarkerElement(; marker, markerpoints, markercolor=color, kwargs...)]
end

function legend_elements(::Type{Lines};
                         color=from_default_theme(:linecolor), kwargs...)
    return [LineElement(; linecolor=color, kwargs...)]
end

function legend_elements(::Type{Contour};
                         color=from_default_theme(:linecolor), kwargs...)
    return [LineElement(; linecolor=color, kwargs...)]
end

function legend_elements(::Type{LinesFill};
                         color=from_default_theme(:linecolor), fillalpha=0.15, kwargs...)
	polycolor = to_color((color, fillalpha))
	return [PolyElement(; polycolor, kwargs...), LineElement(; linecolor=color, kwargs...)]
end

function legend_elements(::Any;
                         color=from_default_theme(:patchcolor), kwargs...)
    return [PolyElement(; polycolor=color, kwargs...)]
end


text:
include("utils.jl")

# starting with tests of simple manifolds
include("euclidean.jl")
include("sphere.jl")
include("rotations.jl")
include("product_manifold.jl")

include("metric_test.jl")


text:
const TO = TrajectoryOptimization

# Get the model
model = Dynamics.car
n,m = model.n, model.m
Q = (1e-2)*Diagonal(I,n)
Qf = 1000.0*Diagonal(I,n)
R = (1e-2)*Diagonal(I,m)
x0 = [0.;0.;0.]
xf = [0.;1.;0.]

obj = TrajectoryOptimization.LQRObjective(Q,R,Qf,xf,N)

# Test Expansion
Q = TO.Expansion{Float64}(n,m)
Q.x .= ones(n)
Q.xx .= ones(n,n)
2*Q
@test Q.x == ones(n)*2
Q/2
@test Q.xx == ones(n,n)

Q2 = copy(Q)
Q2.uu .+= 1
@test Q2.uu == ones(m,m)
@test Q.uu == zeros(m,m)
TO.reset!(Q2)
@test Q2.uu == zeros(m,m)

# Quadratic Cost
Q = Diagonal(1.0I,n)
R = Diagonal(1.0I,m)
H = zeros(m,n)
q = ones(n)
r = zeros(m)
Qf = Diagonal(I,n)*100.
qf = zeros(n)
quadcost = QuadraticCost(Q, R, H, q, r, 0.)
quadcost_f = QuadraticCost(Qf, qf, 0.)
@test quadcost.c == 0
@test quadcost.Q == Q
@test quadcost.R isa Diagonal
@test quadcost_f.c == 0
@test quadcost_f.Q == Qf
@test quadcost_f.q == qf

quadcost = QuadraticCost(Q,R)
@test quadcost.R isa Diagonal
@test quadcost.H == zeros(m,n)

# LQR Cost
xf = [0,1,0]
lqrcost = LQRCost(Q,R,xf)
@test lqrcost.Q == Q
@test lqrcost.q == -Q*xf

cost_term = TO.LQRCostTerminal(Qf, xf)
@test cost_term.Q == Qf
@test cost_term.q == -Qf*xf

# Test costs
x,u = rand(n), rand(m)
dt = rand(1)[1]
@test TO.stage_cost(quadcost,x,u,dt) == 0.5*(x'Q*x + u'R*u)*dt
@test TO.stage_cost(cost_term,x) ≈ 0.5*(x-xf)'Qf*(x-xf)

# Test expansion functions
E = TO.Expansion{Float64}(n,m)
TO.cost_expansion!(E, quadcost, x, u, dt)
@test E.x == Q*x*dt
@test E.xx == Q*dt
@test E.ux == zeros(m,n)*dt
@test u == R*u

TO.reset!(E)
TO.cost_expansion!(E, cost_term, x)
@test E.x ≈ Qf*(x-xf)
@test E.u == zeros(m)
@test E.xx == Qf

grad = PartedVector(zeros(n+m), create_partition((n,m),(:x,:u)))
TO.gradient!(grad, quadcost, x, u, dt)
@test grad.x == dt*Q*x
@test grad.u == dt*R*u

grad = zeros(n)
TO.gradient!(grad, cost_term, x)
@test grad ≈ Qf*(x-xf)

quadcost2 = copy(quadcost)
quadcost2.Q .+= Diagonal(I,n)
@test quadcost2.Q == Diagonal(2I,n)
@test quadcost.Q == Q

# Generic Costs
n,m = 2,1
x,u = rand(n), rand(m)
function mycost(x,u)
    R = Diagonal(0.1I,1)
    Q = 0.1
    return cos(x[1]) + u'R*u + Q*x[2]^2
end
function mycost(xN)
    return cos(xN[1]) + xN[2]^2
end
nlcost = GenericCost(mycost,mycost,n,m)
@test TO.stage_cost(nlcost,x,u) == mycost(x,u)
@test TO.stage_cost(nlcost,x) == mycost(x)

R = Diagonal(0.1I,m)
Q = 0.1
function hess(x,u)
    n,m = length(x),length(u)
    Qexp = Diagonal([-cos(x[1]), 2Q])
    Rexp = 2R
    Hexp = zeros(m,n)
    return Qexp,Rexp,Hexp
end
function hess(x)
    return Diagonal([-cos(x[1]), 2])
end
function gradient(x,u)
    q = [-sin(x[1]), 2Q*x[2]]
    r = 2R*u
    return q,r
end
function gradient(x)
    return [-sin(x[1]), 2*x[2]]
end
nlcost2 = GenericCost(mycost, mycost, gradient, hess, n, m)

E = TO.Expansion{Float64}(n,m)
TO.cost_expansion!(E, nlcost, x, u)
@test E.x == gradient(x,u)[1]
@test E.xx == hess(x,u)[1]
@test E.ux == hess(x,u)[3]
TO.cost_expansion!(E, nlcost2, x, u)
@test E.x == gradient(x,u)[1]
@test E.xx == hess(x,u)[1]
@test E.ux == hess(x,u)[3]

nlcost3 = copy(nlcost)
@test TO.stage_cost(nlcost3,x,u) == TO.stage_cost(nlcost,x,u)

# Minimum time
n = 4
m = 2
N = 10
Q = Diagonal(ones(n))
R = Diagonal(ones(m))
Qf = 0.7*Diagonal(ones(n))
xf = rand(n)

_x = rand(n+1)
_u = rand(m+1)
dt = 0.134

part_mt = (x=1:(n+1), u=(n+1) .+ (1:(m+1)))
part2_mt = (xx=(part_mt.x, part_mt.x), uu=(part_mt.u, part_mt.u), ux=(part_mt.u, part_mt.x), xu=(part_mt.x, part_mt.u))
_cost = LQRCost(Q,R,xf)
_cost_term = TO.LQRCostTerminal(Qf,xf)
_cost_mt = TO.MinTimeCost(_cost,1.6)
_cost_term_mt = TO.MinTimeCost(_cost_term,1.6)

EE = TO.Expansion{Float64}(n+1,m+1)
grad_mt = zeros(n+1+m+1)
grad_term_mt = zeros(n+1)
hess_mt = PartedMatrix(zeros(n+1+m+1,n+1+m+1),part2_mt)
hess_term_mt = zeros(n+1,n+1)
TO.cost_expansion!(EE, _cost_mt, _x, _u, dt)
TO.gradient!(grad_mt,_cost_mt,_x,_u,dt)
TO.gradient!(grad_term_mt,_cost_term_mt,_x)
TO.hessian!(hess_mt,_cost_mt,_x,_u,dt)
TO.hessian!(hess_term_mt,_cost_term_mt,_x)

@test isapprox(EE.x,grad_mt[1:(n+1)])
@test isapprox(EE.u,grad_mt[(n+1) .+ (1:(m+1))])
@test isapprox(EE.xx,hess_mt.xx)
@test isapprox(EE.uu,hess_mt.uu)
@test isapprox(EE.ux,hess_mt.ux)

max_issues_repo_name: MehrdadZareian/GMARA

max_forks_repo_name: MehrdadZareian/GMARA


text:
using YaoBase, YaoArrayRegister

"""
    apply!(register, block)

Apply a block (of quantum circuit) to a quantum register.
"""
function apply!(r::AbstractRegister, b::AbstractBlock)
    _check_size(r, b)
    _apply!(r, b)
end

function _apply!(r::AbstractRegister, b::AbstractBlock)
    _apply_fallback!(r, b)
end

_apply_fallback!(r::AbstractRegister, b::AbstractBlock) =
    throw(NotImplementedError(:_apply_fallback!, (r, b)))

function _apply_fallback!(r::AbstractArrayReg{D,T}, b::AbstractBlock) where {D,T}
    r.state .= mat(T, b) * r.state
    return r
end

"""
    |>(register, blk)

Pipe operator for quantum circuits.

# Example

```julia
julia> arrayreg(bit"0") |> X |> Y
```

!!! warning

    `|>` is equivalent to [`apply!`](@ref), which means it has side effects. You
    need to copy original register, if you do not want to change it in-place.
"""
Base.:(|>)(r::AbstractRegister, blk::AbstractBlock) = apply!(r, blk)

function apply!(r::AbstractRegister, blk::Function)
    if applicable(blk, nactive(r))
        return apply!(r, blk(nactive(r)))
    else
        error(
            "input function is not applicable, it should take a integer as number of current active qubits.",
        )
    end
end

"""
    occupied_locs(x)

Return a tuple of occupied locations of `x`.
"""
occupied_locs(x::AbstractBlock) = (1:nqudits(x)...,)

"""
    subblocks(x)

Returns an iterator of the sub-blocks of a composite block. Default is empty.
"""
subblocks(x::AbstractBlock)

"""
    chsubblocks(composite_block, itr)

Change the sub-blocks of a [`CompositeBlock`](@ref) with given iterator `itr`.
"""
chsubblocks(x::AbstractBlock, itr)

"""
    applymatrix(g::AbstractBlock) -> Matrix

Transform the apply! function of specific block to dense matrix.
"""
applymatrix(T, g::AbstractBlock) =
    linop2dense(T, r -> statevec(apply!(ArrayReg(r), g)), nqudits(g))
applymatrix(g::AbstractBlock) = applymatrix(ComplexF64, g)
# just use BlockMap maybe? No!

print_block(io::IO, blk::AbstractBlock) = print_block(io, MIME("text/plain"), blk)
print_block(blk::AbstractBlock) = print_block(stdout, blk)
print_block(io::IO, ::MIME"text/plain", blk::AbstractBlock) = summary(io, blk)

# return itself by default
Base.copy(x::AbstractBlock) = x

function cache_key end

"""
    mat([T=ComplexF64], blk)

Returns the matrix form of given block.
"""
mat(x::AbstractBlock) = mat(promote_type(ComplexF64, parameters_eltype(x)), x)

mat_matchreg(reg::AbstractRegister, x::AbstractBlock) = mat(x)
mat_matchreg(reg::AbstractArrayReg{D,T}, x::AbstractBlock) where {D,T} = mat(T,x)

Base.Matrix{T}(x::AbstractBlock) where {T} = Matrix(mat(T, x))

# YaoBase interface
YaoBase.nqubits(block::AbstractBlock{2}) = nqudits(block)
YaoBase.nlevel(::Type{<:AbstractBlock{D}}) where {D} = D
YaoBase.nlevel(x::AbstractBlock{D}) where {D} = nlevel(typeof(x))

# properties
for each_property in [:isunitary, :isreflexive, :ishermitian]
    @eval YaoBase.$each_property(x::AbstractBlock) = $each_property(mat(x))
    @eval YaoBase.$each_property(::Type{T}) where {T<:AbstractBlock} =
        $each_property(mat(T))
end

function iscommute_fallback(op1::AbstractBlock{D}, op2::AbstractBlock{D}) where {D}
    _check_block_sizes(op1, op2)
    if length(intersect(occupied_locs(op1), occupied_locs(op2))) == 0
        return true
    else
        return iscommute(mat(op1), mat(op2))
    end
end

YaoBase.iscommute(op1::AbstractBlock{D}, op2::AbstractBlock{D}) where {D} =
    iscommute_fallback(op1, op2)

# parameters
"""
    getiparams(block)

Returns the intrinsic parameters of node `block`, default is an empty tuple.
"""
getiparams(x::AbstractBlock) = ()

"""
    setiparams!([f], block, itr)
    setiparams!([f], block, params...)

Set the parameters of `block`.
When `f` is provided, set parameters of `block` to the value in `collection` mapped by `f`.
`iter` can be an iterator or a symbol, the symbol can be `:zero`, `:random`.
"""
function setiparams! end

"""
    setiparams([f], block, itr)
    setiparams([f], block, params...)

Set the parameters of `block`, the non-inplace version.
When `f` is provided, set parameters of `block` to the value in `collection` mapped by `f`.
`iter` can be an iterator or a symbol, the symbol can be `:zero`, `:random`.
"""
function setiparams end

for F in [:setiparams!, :setiparams]
    @eval begin
        $F(x::AbstractBlock, args...) =
            niparams(x) == length(args) == 0 ? x :
            throw(NotImplementedError($(QuoteNode(F)), (x, args...)))

        $F(x::AbstractBlock, it::Union{Tuple,AbstractArray,Base.Generator}) = $F(x, it...)
        $F(x::AbstractBlock, a::Number, xs::Number...) =
            error("setparams!(x, θ...) is not implemented")
        $F(x::AbstractBlock, it::Symbol) = $F(x, render_params(x, it))

        $F(f::Function, x::AbstractBlock, it) =
            $F(x, map(x -> f(x...), zip(getiparams(x), it)))
        $F(f::Nothing, x::AbstractBlock, it) = $F(x, it)
        $F(f::Function, x::AbstractBlock, it::Symbol) = $F(f, x, render_params(x, it))
    end
end

"""
    parameters(block)

Returns all the parameters contained in block tree with given root `block`.
"""
parameters(x::AbstractBlock) = parameters!(parameters_eltype(x)[], x)

"""
    parameters!(out, block)

Append all the parameters contained in block tree with given root `block` to
`out`.
"""
function parameters!(out, x::AbstractBlock)
    append!(out, getiparams(x))
    for blk in subblocks(x)
        parameters!(out, blk)
    end
    return out
end

# = prewalk(blk->append!(out, getiparams(blk)), x)

niparams(x::AbstractBlock) = length(getiparams(x))

function nparameters(x::AbstractBlock)
    count = niparams(x)
    for each in subblocks(x)
        count += nparameters(each)
    end
    return count
end

"""
    iparams_eltype(block)

Return the element type of [`getiparams`](@ref).
"""
iparams_eltype(x::AbstractBlock) = eltype(getiparams(x))

"""
    parameters_eltype(x)

Return the element type of [`parameters`](@ref).
"""
function parameters_eltype(x::AbstractBlock)
    T = iparams_eltype(x)
    for each in subblocks(x)
        T = promote_type(T, parameters_eltype(each))
    end
    return T
end

mutable struct Dispatcher{VT}
    params::VT
    loc::Int
end

Dispatcher(params) = Dispatcher(params, 0)

function consume!(d::Dispatcher, n::Int)
    d.loc += n
    d.params[d.loc-n+1:d.loc]
end

function consume!(d::Dispatcher{<:Symbol}, n::Int)
    d.loc += n
    d.params
end

function consume!(d::Dispatcher{<:Number}, n::Int)
    if n == 0
        return ()
    elseif n == 1
        d.loc += n
        return d.params
    else
        error("do not have enough parameters to consume, expect 0, 1, got $n")
    end
end

function dispatch!(f::Union{Function,Nothing}, x::AbstractBlock, it::Dispatcher)
    setiparams!(f, x, consume!(it, niparams(x)))
    for each in subblocks(x)
        dispatch!(f, each, it)
    end
    return x
end

"""
    dispatch!(x::AbstractBlock, collection)

Dispatch parameters in collection to block tree `x`.

!!! note

    it will try to dispatch the parameters in collection first.
"""
function dispatch!(f::Union{Function,Nothing}, x::AbstractBlock, it)
    dp = Dispatcher(it)
    res = dispatch!(f, x, dp)
    @assert (it isa Symbol || length(it) == dp.loc) "expect $(dp.loc) parameters, got $(length(it))"
    return res
end

dispatch!(x::AbstractBlock, it) = dispatch!(nothing, x, it)

"""
    popdispatch!(f, block, list)

Pop the first [`nparameters`](@ref) parameters of list, map them with a function
`f`, then dispatch them to the block tree `block`. See also [`dispatch!`](@ref).
"""
function popdispatch!(f::Function, x::AbstractBlock, list::Vector)
    setiparams!(f, x, (popfirst!(list) for k = 1:niparams(x))...)
    for each in subblocks(x)
        popdispatch!(f, each, list)
    end
    return x
end

"""
    popdispatch!(block, list)

Pop the first [`nparameters`](@ref) parameters of list, then dispatch them to
the block tree `block`. See also [`dispatch!`](@ref).
"""
function popdispatch!(x::AbstractBlock, list::Vector)
    setiparams!(x, (popfirst!(list) for k = 1:niparams(x))...)
    for each in subblocks(x)
        popdispatch!(each, list)
    end
    return x
end

render_params(r::AbstractBlock, params) = params
render_params(r::AbstractBlock, params::Symbol) = render_params(r, Val(params))
render_params(r::AbstractBlock, ::Val{:random}) = (rand() for i = 1:niparams(r))
render_params(r::AbstractBlock, ::Val{:zero}) =
    (zero(iparams_eltype(r)) for i = 1:niparams(r))

"""
    cache_type(::Type) -> DataType

Return the element type that a [`CacheFragment`](@ref)
will use.
"""
cache_type(::Type{<:AbstractBlock}) = Any

"""
    cache_key(block)

Returns the key that identify the matrix cache of this block. By default, we
use the returns of [`parameters`](@ref) as its key.
"""
cache_key(x::AbstractBlock)

function _check_size(r::AbstractRegister{D}, pb::AbstractBlock{D}) where D
    nqudits(pb) == nactive(r) ||
        throw(QubitMismatchError("register size $(nactive(r)) mismatch with block size $(nqudits(pb))"))
    nlevel(r) == nlevel(pb) ||
        throw(QubitMismatchError("register number of level $(nlevel(r)) mismatch with block number of level $(nlevel(pb))"))
end

_check_block_sizes(blocks::AbstractBlock{D}...) where D = _check_block_sizes(blocks)
function _check_block_sizes(blocks::Union{AbstractVector{<:AbstractBlock{D}}, NTuple{N, AbstractBlock{D}} where N}, n=nothing) where D
    if n === nothing
        if length(blocks) == 0
            error("number of blocks can not be 0!")
        else
            n = nqudits(first(blocks))
        end
    end
    for b in blocks
        if nqudits(b) != n
            throw(QubitMismatchError("number of qudits `$(nqudits(b))` does not match the expected value `$(n)`."))
        end
    end
    return n
end

"""
    parameters_range(block)

Return the range of real parameters present in `block`.

!!! note
    It may not be the case that `length(parameters_range(block)) == nparameters(block)`.

# Example

```jldoctest; setup=:(using YaoBlocks)
julia> YaoBlocks.parameters_range(RotationGate(X, 0.1))
1-element Vector{Tuple{Float64, Float64}}:
 (0.0, 6.283185307179586)
```
"""
function parameters_range(block::AbstractBlock)
    T = parameters_eltype(block)
    out = Tuple{T,T}[]
    parameters_range!(out, block)
    out
end

function parameters_range!(out::Vector{Tuple{T,T}}, block::AbstractBlock) where {T}
    np = niparams(block)
    if  np!= 0
        ranges = iparams_range(block)
        @assert length(ranges) == np
        for i=1:np
            push!(out, ranges[i])
        end
    end
    for subblock in subblocks(block)
        parameters_range!(out, subblock)
    end
end

# non-inplace versions
"""
    apply(register, block)

The non-inplace version of applying a block (of quantum circuit) to a quantum register.
Check `apply!` for the faster inplace version.
"""
apply(r::AbstractRegister, b) = apply!(copy(r), b)

function generic_dispatch!(f::Union{Function,Nothing}, x::AbstractBlock, it::Dispatcher)
    x = setiparams(f, x, consume!(it, niparams(x)))
    chsubblocks(x, map(subblocks(x)) do each
        generic_dispatch!(f, each, it)
    end)
end

"""
    dispatch(x::AbstractBlock, collection)

Dispatch parameters in collection to block tree `x`, the generic non-inplace version.

!!! note

    it will try to dispatch the parameters in collection first.
"""
function dispatch(f::Union{Function,Nothing}, x::AbstractBlock, it)
    dp = Dispatcher(it)
    res = generic_dispatch!(f, x, dp)
    @assert (it isa Symbol || length(it) == dp.loc) "expect $(dp.loc) parameters, got $(length(it))"
    return res
end

dispatch(x::AbstractBlock, it) = dispatch(nothing, x, it)


text:
function housenumberssum(inputarray)
    accu = 0
    for x in inputarray
        if x == 0 return accu end
        accu += x
    end
end

max_stars_repo_name: UlrichBerntien/Codewars-Katas

max_issues_repo_name: UlrichBerntien/Codewars-Katas

max_forks_repo_name: UlrichBerntien/Codewars-Katas


text:
L = 1:1000; d = 1:1000;
L2d = Array{Float64}(undef, length(L), length(d));
c(L,d) = (10*L/(2*pi)*tanh(2*pi*d/L))^0.5;

Lmatrix = [fill(convert(Float64,i),(1,length(d))) for i in 1:length(d)]


text:
#===========================================
Transport matrices
===========================================#
T_D(p) = T_Circulation
S₀ = PFDO(ones(nb), DIV, Iabove)
S′ = PFDO(ztop, DIV, Iabove)
function T_POP(p)
    @unpack w₀, w′ = p
    return w₀ * S₀ + w′ * S′
end
T_all = (T_D, T_D, T_POP)

#===========================================
Sources minus sinks
===========================================#
# Geological Restoring
function geores(x, p)
    @unpack τgeo, xgeo = p
    return @. (xgeo - x) / τgeo
end
# Uptake of phosphate (DIP)
relu(x) = (x .≥ 0) .* x
function uptake(DIP, p)
    @unpack τDIP, k = p
    DIP⁺ = relu(DIP)
    return @. 1 / τDIP * DIP⁺^2 / (DIP⁺ + k) * (ztop == 0)
end
# Remineralization DOP into DIP
function remineralization(DOP, p)
    @unpack τDOP = p
    return @. DOP / τDOP
end
# Dissolution of POP into DOP
function dissolution(POP, p)
    @unpack τPOP = p
    return @. POP / τPOP
end
# Add them up into sms functions (Sources Minus Sinks)
function sms_DIP(DIP, DOP, POP, p)
    return -uptake(DIP, p) + remineralization(DOP, p) + geores(DIP, p)
end
function sms_DOP(DIP, DOP, POP, p)
    @unpack σ = p
    return σ * uptake(DIP, p) - remineralization(DOP, p) + dissolution(POP, p)
end
function sms_POP(DIP, DOP, POP, p)
    @unpack σ = p
    return (1 - σ) * uptake(DIP, p) - dissolution(POP, p)
end
sms_all = (sms_DIP, sms_DOP, sms_POP) # bundles all the source-sink functions in a tuple

#===========================================
In-place sources minus sinks
===========================================#

function G_DIP!(dx, DIP, DOP, POP, p)
    @unpack τgeo, xgeo, τDIP, k, τDOP = p
    @. dx = (xgeo - DIP) / τgeo - (DIP ≥ 0) / τDIP * DIP^2 / (DIP + k) * (ztop == 0) + DOP / τDOP
end
function G_DOP!(dx, DIP, DOP, POP, p)
    @unpack τgeo, xgeo, τDIP, k, τDOP, τPOP, σ = p
    @. dx = σ * (DIP ≥ 0) / τDIP * DIP^2 / (DIP + k) * (ztop == 0) - DOP / τDOP + POP / τPOP
end
function G_POP!(dx, DIP, DOP, POP, p)
    @unpack τgeo, xgeo, τDIP, k, τDOP, τPOP, σ = p
    @. dx = (1 - σ) * (DIP ≥ 0) / τDIP * DIP^2 / (DIP + k) * (ztop == 0) - POP / τPOP
end
Gs = (G_DIP!, G_DOP!, G_POP!)

#===========================================
AIBECS F and ∇ₓF
===========================================#
fun = AIBECSFunction(T_all, sms_all, nb)
fun! = AIBECSFunction((T_D, T_POP), Gs, nb, 3, [1, 1, 2])
F, ∇ₓF = F_and_∇ₓF(T_all, sms_all, nb) # deprecated
F!, ∇ₓF! = F_and_∇ₓF((T_D, T_POP), Gs, nb, 3, [1, 1, 2]) # deprecated

#===========================================
AIBECS split operators for CNLF
===========================================#

NL_DIP(DIP, DOP, POP, p) = -uptake(DIP, p) + geores(0 * DIP, p)
L_DIP(DIP, DOP, POP, p) = remineralization(DOP, p) + geores(DIP, p) - geores(0 * DIP, p)

function NL_DOP(DIP, DOP, POP, p)
    @unpack σ = p
    return σ * uptake(DIP, p)
end
L_DOP(DIP, DOP, POP, p) = -remineralization(DOP, p) + dissolution(POP, p)

function NL_POP(DIP, DOP, POP, p)
    @unpack σ = p
    return (1 - σ) * uptake(DIP, p)
end
L_POP(DIP, DOP, POP, p) = -dissolution(POP, p)

NL_all = (NL_DIP, NL_DOP, NL_POP)
L_all = (L_DIP, L_DOP, L_POP)

F2, L, NL, ∇ₓF2, ∇ₓL, ∇ₓNL, T = split_state_function_and_Jacobian(T_all, L_all, NL_all, nb)

#===========================================
Tests
===========================================#

@testset "Biogeochemical functions" begin
    @testset "Particle Flux Divergence is conservative" begin
        @test norm(v) / norm(T_POP(p)' * v) > ustrip(upreferred(1u"Myr"))
    end
    nt = length(T_all)
    n = nt * nb
    @unpack xgeo = p
    x = xgeo * ones(n) .* exp.(cos.(collect(1:n)) / 10)
    @testset "State function" begin
        F₀ = fun(x, p)
        v_all = repeat(v, nt)
        @test F₀ isa Vector{Float64}
        @test size(F₀) == (n,)
        @test norm(v_all .* x) / norm(v_all' * F₀) > ustrip(upreferred(1u"Myr"))
    end
    @testset "Jacobian of the state function" begin
        ∇ₓF₀ = fun.jac(x, p)
        @test ∇ₓF₀ isa SparseMatrixCSC
        @test size(∇ₓF₀) == (n, n)
    end
end

@testset "In-place functions" begin
    nt = length(T_all)
    n = nt * nb
    @unpack xgeo = p
    x = xgeo * ones(n) .* exp.(cos.(collect(1:n)) / 10)
    dx = similar(x)
    @testset "in-place F ≈ out-of-place F" begin
        @test fun(x, p) ≈ fun!(dx, x, p) rtol = 1e-10
        @test fun.jac(x, p) ≈ fun!.jac(x, p) rtol = 1e-10
    end
end

@testset "Linear/nonlinear split functions" begin
    nt = length(T_all)
    n = nt * nb
    @unpack xgeo = p
    x = xgeo * ones(n) .* exp.(cos.(collect(1:n)) / 10)
    @testset "split F ≈ F" begin
        @test F2(x, p) ≈ fun(x, p) rtol = 1e-10
    end
    @testset "split F ≈ L + NL" begin
        @test F2(x, p) ≈ L(x, p) + NL(x, p) - T(p) * x rtol = 1e-10
    end
    @testset "split ∇ₓF ≈ ∇ₓF" begin
        @test fun.jac(x, p) ≈ ∇ₓF2(x, p) rtol = 1e-10
    end
    @testset "split ∇ₓF ≈ T + ∇ₓNL" begin
        @test ∇ₓF2(x, p) ≈ ∇ₓL(p) + ∇ₓNL(x, p) - T(p) rtol = 1e-10
    end
end


text:
module TestVariantVisualization

using VariantVisualization
using Test
using DataFrames
using GeneticVariation
#using Base.Test #version
using DelimitedFiles

all_tests = [
    ("new_vcf_utils.jl",   "           Testing: New VCF Utils"),
#    ("plot_utils.jl",     "       Testing: Plot Utils")
    ]

println("Running tests:")

for (t, test_string) in all_tests
    println("-----------------------------------------")
    println("-----------------------------------------")
    println(test_string)
    println("-----------------------------------------")
    println("-----------------------------------------")
    include(t)
end

end


text:
using ElectromagneticFields
using Test

println()
include("test_analytic.jl")


text:
using Plots, SFFM

# include("../../src/SFFM.jl")

cme_9 = SFFM.MakeME(SFFM.CMEParams[9])

f = SFFM.pdf(cme_9)
F(x) = 1 - SFFM.cdf(cme_9)(x)

x = 0:0.05:1.5
plot(x,f.(x), label = "α exp(Sz) s")
plot!(x,f.(x.+0.3)./F.(0.3), label = "α exp(S(0.3+z)) s/α exp(S 0.3) e")
plot!(x,f.(x.+0.6)./F.(0.6), label = "α exp(S(0.6+z)) s/α exp(S 0.6) e")
plot!(
    xlabel = "z",
    ylabel = "Density",
    legend = :outertop,
)

savefig("examples/meNumerics/ME_residual_life_density.pdf")

max_stars_repo_name: angus-lewis/SFFM

max_issues_repo_name: angus-lewis/SFFM

max_forks_repo_name: angus-lewis/SFFM


text:
# ---
# title: 1535. Find the Winner of an Array Game
# id: problem1535
# author: Tian Jun
# date: 2020-10-31
# difficulty: Medium
# categories: Array
# link: <https://leetcode.com/problems/find-the-winner-of-an-array-game/description/>
# hidden: true
# ---
#
# Given an integer array `arr` of **distinct** integers and an integer `k`.
#
# A game will be played between the first two elements of the array (i.e.
# `arr[0]` and `arr[1]`). In each round of the game, we compare `arr[0]` with
# `arr[1]`, the larger integer wins and remains at position `0` and the smaller
# integer moves to the end of the array. The game ends when an integer wins `k`
# consecutive rounds.
#
# Return _the integer which will win the game_.
#
# It is **guaranteed** that there will be a winner of the game.
#
#
#
# **Example 1:**
#
#
#
#     Input: arr = [2,1,3,5,4,6,7], k = 2
#     Output: 5
#     Explanation: Let's see the rounds of the game:
#     Round |       arr       | winner | win_count
#       1   | [2,1,3,5,4,6,7] | 2      | 1
#       2   | [2,3,5,4,6,7,1] | 3      | 1
#       3   | [3,5,4,6,7,1,2] | 5      | 1
#       4   | [5,4,6,7,1,2,3] | 5      | 2
#     So we can see that 4 rounds will be played and 5 is the winner because it wins 2 consecutive games.
#
#
# **Example 2:**
#
#
#
#     Input: arr = [3,2,1], k = 10
#     Output: 3
#     Explanation: 3 will win the first 10 rounds consecutively.
#
#
# **Example 3:**
#
#
#
#     Input: arr = [1,9,8,2,3,7,6,4,5], k = 7
#     Output: 9
#
#
# **Example 4:**
#
#
#
#     Input: arr = [1,11,22,33,44,55,66,77,88,99], k = 1000000000
#     Output: 99
#
#
#
#
# **Constraints:**
#
#   * `2 <= arr.length <= 10^5`
#   * `1 <= arr[i] <= 10^6`
#   * `arr` contains **distinct**  integers.
#   * `1 <= k <= 10^9`
#
#
## @lc code=start
using LeetCode

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


text:
using DataTypesBasic
using Documenter

makedocs(;
    modules=[DataTypesBasic],
    authors="Stephan Sahm <stephan.sahm@gmx.de> and contributors",
    repo="https://github.com/JuliaFunctional/DataTypesBasic.jl/blob/{commit}{path}#L{line}",
    sitename="DataTypesBasic.jl",
    format=Documenter.HTML(;
        prettyurls=get(ENV, "CI", "false") == "true",
        canonical="https://JuliaFunctional.github.io/DataTypesBasic.jl",
        assets=String[],
    ),
    pages=[
        "Home" => "index.md",
        "Manual" => "manual.md",
        "Library" => "library.md",
    ],
)

deploydocs(;
    repo="github.com/JuliaFunctional/DataTypesBasic.jl",
    devbranch="main",
)