search

SciML / DiffEqBase.jl

The lightweight Base library for shared types and functionality for defining differential equation and scientific machine learning (SciML) problems
Other
313 stars 112 forks source link

ODE_DEFAULT_NORM fails for ArrayPartition containing empty arrays #664

Closed mfsch closed 3 years ago

mfsch commented 3 years ago

When one of the arrays making up an array partition is empty, calls to ODE_DEFAULT_NORM fail:

julia> u = ArrayPartition(zeros(1), zeros(0))
([0.0], Float64[])

julia> DiffEqBase.ODE_DEFAULT_NORM(u, nothing)
ERROR: ArgumentError: reducing over an empty collection is not allowed

I’m not sure whether #619 & #620 have introduced this behavior or if it existed before (I’m pretty sure it used to work in some much older version), but it appears that sum(UNITLESS_ABS, x) is called on each array of the partition and fails for empty arrays. A possible solution would be to include an init=zero(eltype(x)), but that only works for Julia 1.6 so it’s not ideal. After a bit of digging, another option appears to be adding a specialization for mapreduce_empty:

Base.mapreduce_empty(::typeof(UNITLESS_ABS), op, T) = abs2(Base.reduce_empty(op, T))

The context in which I’m running into this issue is that I’m solving PDEs with variables on a staggered grid that’s distributed over MPI, so one of the variables ends up having no values on one of the MPI ranks.

ChrisRackauckas commented 3 years ago

I think specializing mapreduce_empty would make sense here.