Add functionality to keep singleton dimensions for reduction

In Julia, functions like `mean` keep the singleton dimension. In ClimaAnalysis, we don't follow this because it is ambiguous what the value should be. However, there should be an option to keep the singleton dimension. For example, the keyword `keepdim` can be added to the function. For the value of the dimension, it could default to the first value of the dimension.

Start by locating ClimaAnalysis reduction functions and their existing dimension-handling tests. Review how Julia's reduction functions preserve singleton dimensions, then define the keepdim option and its default dimension value. Done means supported reductions preserve the singleton dimension when requested and existing behavior remains unchanged otherwise.

Written by the indexing model from the issue text.