Remove output dimension type parameter and add output_dim to API

[devmotion]

The PR removes the type parameter. I tested with my downstream package and tests finally pass again with this change.

The reason is that

• in general, the output dimension cannot be inferred from u
• this implies that the constructor of basically all interpolation methods can't be inferred
• I think it's questionable how useful the information is anyway given that output dimension 3 can refer to u::Vector{<:Vector} where each element is of length 3 or e.g. a u::Matrix with 3 rows

On the master branch:

julia> using DataInterpolations, Test

julia> @inferred(ConstantInterpolation([1;; 2;; 3], [0.1, 0.2, 0.3]))
ERROR: return type ConstantInterpolation{Matrix{Int64}, Vector{Float64}, Vector{Union{}}, Int64, (1,)} does not match inferred return type ConstantInterpolation{Matrix{Int64}, Vector{Float64}, Vector{T}, Int64} where T
Stacktrace:
 [1] error(s::String)
   @ Base ./error.jl:35
 [2] top-level scope
   @ REPL[13]:1

julia> @inferred(ConstantInterpolation(["A", "B", "C"], [0.1, 0.2, 0.3]))
ERROR: return type ConstantInterpolation{Vector{String}, Vector{Float64}, Nothing, String, (1,)} does not match inferred return type ConstantInterpolation{Vector{String}, Vector{Float64}, Nothing, String}
Stacktrace:
 [1] error(s::String)
   @ Base ./error.jl:35
 [2] top-level scope
   @ REPL[14]:1

julia> @inferred(ConstantInterpolation([[1], [2], [3]], [0.1, 0.2, 0.3]))
ERROR: return type ConstantInterpolation{Vector{Vector{Int64}}, Vector{Float64}, Nothing, Vector{Int64}, (1,)} does not match inferred return type ConstantInterpolation{Vector{Vector{Int64}}, Vector{Float64}, Nothing, Vector{Int64}}
Stacktrace:
 [1] error(s::String)
   @ Base ./error.jl:35
 [2] top-level scope
   @ REPL[15]:1

This PR:

julia> using DataInterpolations, Test

julia> @inferred(ConstantInterpolation([1;; 2;; 3], [0.1, 0.2, 0.3]))
ConstantInterpolation with 3 points, in left direction
┌──────┬─────┐
│ time │  u1 │
├──────┼─────┤
│  0.1 │ 1.0 │
│  0.2 │ 2.0 │
│  0.3 │ 3.0 │
└──────┴─────┘


julia> @inferred(ConstantInterpolation(["A", "B", "C"], [0.1, 0.2, 0.3]))
ConstantInterpolation with 3 points, in left direction
┌──────┬───┐
│ time │ u │
├──────┼───┤
│  0.1 │ A │
│  0.2 │ B │
│  0.3 │ C │
└──────┴───┘


julia> @inferred(ConstantInterpolation([[1], [2], [3]], [0.1, 0.2, 0.3]))
ConstantInterpolation with 3 points, in left direction
┌──────┬─────┐
│ time │  u1 │
├──────┼─────┤
│  0.1 │ 1.0 │
│  0.2 │ 2.0 │
│  0.3 │ 3.0 │
└──────┴─────┘

---

Edit: I added an exported output_dim function as discussed in #396 (comment). This could be used by SparseConnectivityTracer if desired/useful which currently seems to be the only downstream user of the type parameter but only makes us of it in a non-static way.

[SouthEndMusic]

@ChrisRackauckas isn't this another argument in favor of the SimpleInterpolations.jl you suggested?

I guess this feature has some use somewhere otherwise it wouldn't have been added

[devmotion]

IMO it's just generally bad and not useful: It can't be inferred and doesn't tell you anything about the structure of u (and hence the output type).

And the type inference problem implies that generally you pay runtime dispatch penalties in a function if you construct the interpolation and use it for computations.

I guess this feature has some use somewhere otherwise it wouldn't have been added

Downstream packages could still use the same runtime logic for deriving the output dimension (could even be made available as part of the API) if they want that information, without breaking type inference for every user of DataInterpolations.

[sathvikbhagavan]

I had added it as a feature request (#328) but if it breaks inference, I guess its not worth it.

[devmotion]

It's used https://github.com/adrhill/SparseConnectivityTracer.jl/blob/622cfb729832d127f2d0f9f94110391b4d439b6b/ext/SparseConnectivityTracerDataInterpolationsExt.jl#L54 and https://github.com/adrhill/SparseConnectivityTracer.jl/blob/622cfb729832d127f2d0f9f94110391b4d439b6b/ext/SparseConnectivityTracerDataInterpolationsExt.jl#L64 but these cases do not need a static type information but can just extract the dimension (IMO still questionable how useful in general given the different types of u).

[ChrisRackauckas]

Why is it not inferable? This seems trivial to infer?

[devmotion]

No, it can only be inferred from static arrays. It queries length(first(u)) and size(u, ndims(u)), depending on the type of u. And the resulting number doesn't even let you distinguish between the two cases, which return completely different output types when evaluating multiple time points.

[ChrisRackauckas]

Wait it's not a dimension size computable from ndims(u)? I guess I haven't looked at it in awhile but it should be as simple as that.

[devmotion]

No, it's not. It refers to the actual size of the output/input data.

[ChrisRackauckas]

Oh, that can't be type information. Wait why would DI even want or need that? That seems silly.

[devmotion]

🤷

[SouthEndMusic]

From the original issue by @adrhill:

As a developer writing a package extension for DataInterpolations.jl (adrhill/SparseConnectivityTracer.jl#178), it would be useful if all AbstractInterpolations indicated their output dimensionality, e.g.

So it might be useful to use get_output_dim in the DI extension directly, so that DI doesn't have to know explicitly how different types of u are sliced by DataInterpolations.jl? Or should they just infer it from an evaluation?

[adrhill]

Wait it's not a dimension size computable from ndims(u)? I guess I haven't looked at it in awhile but it should be as simple as that.

I thought this as well. At least that was my intention in opening #328 – Having the interpolations mirror array types. I have no problem with a revert/fix of whatever N is currently doing.

[ChrisRackauckas]

Let's fix it to ndims. This seems to just be a mistake.

[devmotion]

So, you'd want

output_dims(::AbstractVector) = 0 # each value is a scalar
output_dims(::AbstractVector{<:AbstractArray{<:Any,N}}) where {N} = N # each value is an array but values are not stacked
output_dims(::AbstractArray{<:Any,N}) where {N} = N - 1 # each value is an array but multiple values are stacked

Note that you still wouldn't be able to draw any conclusion about the structure of the output based on N.