(Abstract) reductions don't play well with other operations

[tomchor]

Currently I can't operate on a reduction using AbstractOperations:

julia> using Oceananigans;

julia> grid = RectilinearGrid(size=(4,4,4,), extent=(1,1,1));

julia> model = NonhydrostaticModel(grid = grid);

julia> Average(model.velocities.u)/2
ERROR: MethodError: no method matching /(::Reduction{typeof(Statistics.mean!), Field{Face, Center, Center, Nothing, RectilinearGrid{Float64, Periodic, Periodic, Bounded, Float64, Float64, Float64, OffsetArrays.OffsetVector{Float64, StepRangeLen{Float64, Base.TwicePrecision{Float64}, Base.TwicePrecision{Float64}}}, OffsetArrays.OffsetVector{Float64, StepRangeLen{Float64, Base.TwicePrecision{Float64}, Base.TwicePrecision{Float64}}}, OffsetArrays.OffsetVector{Float64, StepRangeLen{Float64, Base.TwicePrecision{Float64}, Base.TwicePrecision{Float64}}}, CPU}, Tuple{Colon, Colon, Colon}, OffsetArrays.OffsetArray{Float64, 3, Array{Float64, 3}}, Float64, FieldBoundaryConditions{BoundaryCondition{Oceananigans.BoundaryConditions.Periodic, Nothing}, BoundaryCondition{Oceananigans.BoundaryConditions.Periodic, Nothing}, BoundaryCondition{Oceananigans.BoundaryConditions.Periodic, Nothing}, BoundaryCondition{Oceananigans.BoundaryConditions.Periodic, Nothing}, BoundaryCondition{Oceananigans.BoundaryConditions.Flux, Nothing}, BoundaryCondition{Oceananigans.BoundaryConditions.Flux, Nothing}, BoundaryCondition{Oceananigans.BoundaryConditions.Flux, Nothing}}, Nothing, Oceananigans.Fields.FieldBoundaryBuffers{Nothing, Nothing, Nothing, Nothing}}, Tuple{Int64, Int64, Int64}}, ::Int64)
Closest candidates are:
  /(::Any, ::Any, ::Any, ::Oceananigans.Grids.AbstractGrid, ::Any, ::Any, ::Number, ::Oceananigans.AbstractOperations.BinaryOperation) at /glade/work/tomasc/.julia/packages/Oceananigans/Ey1oO/src/AbstractOperations/binary_operations.jl:88
  /(::Any, ::Any, ::Any, ::Oceananigans.Grids.AbstractGrid, ::Any, ::Any, ::Oceananigans.AbstractOperations.BinaryOperation, ::Number) at /glade/work/tomasc/.julia/packages/Oceananigans/Ey1oO/src/AbstractOperations/binary_operations.jl:83
  /(::Any, ::Any, ::Any, ::Oceananigans.Grids.AbstractGrid, ::Any, ::Any, ::Oceananigans.AbstractOperations.BinaryOperation, ::Oceananigans.AbstractOperations.BinaryOperation) at /glade/work/tomasc/.julia/packages/Oceananigans/Ey1oO/src/AbstractOperations/binary_operations.jl:67
  ...
Stacktrace:
 [1] top-level scope
   @ REPL[11]:1
 [2] top-level scope
   @ /glade/work/tomasc/.julia/packages/CUDA/DfvRa/src/initialization.jl:52

This is the case for Average and Integral abstract operators (but not maximum and minimum) and any other operation that I tried.

Is this expected?

[glwagner]

Right, we cannot form operations with reductions. To operate on reductions, we have to wrap them in a Field, so that they can be computed, stored, and then used via the stored result. You can write

U = Field(Average(model.velocities.u))
U / 2

[glwagner]

We could to suggest wrapping the Reduction in a Field. I'm not 100% sure the best way to inject that kind of hint; we don't own operations like / so I think to do that we would have to define them (and then throw an error). Hmm.

[navidcy]

Would this work?

julia> Average(model.velocities.u / 2)

[tomchor]

Would this work?

julia> Average(model.velocities.u / 2)

Yup, this works. And that's what I ended up doing. I just thought I'd post it here since I wasn't sure this was expected behavior.

Should I close the issue?

[glwagner]

There's no way around it since we certainly cannot compute reductions on the fly, we have to precompute them. But maybe there is a part of the docs or doc string we can improve? Not sure if this was expected because of some part of the docs

[glwagner]

This behavior is by design so closing