"""Create a `storage` array for [`propagate`](@ref).

```julia
storage = init_storage(state, tlist)
```

creates a storage array suitable for storing a `state` for each point in
`tlist`.

```julia
storage = init_storage(state, tlist, observables))
```

creates a storage array suitable for the data generated by the `observables`
applied to `state`, see [`map_observables`](@ref), for each point in `tlist`.

```julia
storage = init_storage(data, nt))
```

creates a storage arrays suitable for storing `data` nt times, where
`nt=length(tlist)`. By default, this will be a vector of `typeof(data)` and
length `nt`, or a `n × nt` Matrix with the same `eltype` as `data` if `data` is
a Vector of length `n`.
"""
function init_storage(state, tlist::AbstractVector)
    nt = length(tlist)
    return init_storage(state, nt)
end

function init_storage(state, tlist, observables)
    data = map_observables(observables, state)
    nt = length(tlist)
    return init_storage(data, nt)
end

init_storage(data, nt::Integer) = Vector{typeof(data)}(undef, nt)

init_storage(data::Vector, nt::Integer) = Matrix{eltype(data)}(undef, length(data), nt)


"""Obtain "observable" data from `state`.

```julia
data = map_observables(observables, state)
```

calculates the data for a tuple of `observables` applied to `state`.
For a single observable (tuple of length 1), simply return the result of
[`map_observable`](@ref).

For multiple observables, return the tuple resulting from applying
[`map_observable`](@ref) for each observable. If the tuple is "uniform" (all
elements are of the same type, e.g. if each observable calculates the
expectation value of a Hermitian operator), it is converted to a Vector. This
allows for compact storage in a storage array, see [`init_storage`](@ref).
"""
function map_observables(observables, state)
    if length(observables) == 1
        return map_observable(observables[1], state)
    else
        val_tuple = Tuple(map_observable(O, state) for O in observables)
        uniform_type = typeof(val_tuple[1])
        is_uniform = all(typeof(v) == uniform_type for v in val_tuple[2:end])
        if is_uniform
            return collect(val_tuple)  # convert to Vector
        else
            return val_tuple
        end
    end
end


"""Apply a single `observable` to `state`.

```julia
data = map_observable(observable, state)
```

By default, `observable` is assumed to be callable, and the above is equivalent
to `data = observable(state)`.

If `observable` is a matrix and `state` is a vector evaluate the expectation
value of the observable as `dot(state, observable, state)`.

"""
function map_observable(observable, state)
    return observable(state)
end

function map_observable(observable::AbstractMatrix, state::AbstractVector)
    return dot(state, observable, state)
end


"""Place data into `storage` for time slot `i`.

```julia
    write_to_storage!(storage, i, state, observables)
```

For a `storage` array created by [`init_storage`](@ref), store the data obtains
from [`map_observables`](@ref) into the `storage` for time slot `i`. This
delegates to the more general

```julia
    write_to_storage!(storage, i, data)
```

Conceptually, this corresponds roughly to `storage[i] = data`, but `storage`
may have its own idea on how to store data for a specific time slot. For
example, with the default [`init_storage`](@ref) Vector data will be stored in
a matrix, and `write_to_storage!` will in this case write data to the i'th
column of the matrix.

For a given type of `storage` and `data`, it is the developer's responsibility
that `init_storage` and `write_to_storage!` are compatible.
"""
function write_to_storage!(storage, i::Integer, state, observables)
    data = map_observables(observables, state)
    write_to_storage!(storage, i, data)
end

function write_to_storage!(storage::AbstractVector, i::Integer, data)
    storage[i] = data
end

function write_to_storage!(storage::Matrix{T}, i::Integer, data::Vector{T}) where T
    storage[:, i] .= data
end


"""Obtain data from storage

```julia
get_from_storage!(state, storage, i)
```

extracts data from the `storage` for the i'th time slot. Invese of
[`write_to_storage!`](@ref)
"""
get_from_storage!(state, storage::AbstractVector, i) = copyto!(state, storage[i])
get_from_storage!(state, storage::Matrix, i) = copyto!(state, storage[:, i])
