Should lower precision dtypes be inferred in linspace, etc ?

This is essentially the same problem reported in #8597, but a part that I feel was not fully addressed when it was recently closed by a doc improvement. Essentially, linspace, geomspace, arange, r_, etc. will never infer a dtype with less precision than float64 or complex128. It is unclear whether this is intended behavior.

The docs for linspace, arange, logspace, geomspace all have something along the lines of "If dtype is not given, the data type is inferred from [other arguments]." It clarifies that ints will be promoted to floats, but not that float32 will be upgraded to float64s, while float128 will be left alone (similarly with complex64/128/256).

import numpy as np

print([ # same for geomspace, logspace
    np.linspace(np.float16(0.1), np.float16(0.5), 3).dtype,
    np.linspace(np.float32(0.1), np.float32(0.5), 3).dtype,
    np.linspace(np.float64(0.1), np.float64(0.5), 3).dtype,
    np.linspace(np.float128(0.1), np.float128(0.5), 3).dtype
    ])

print([
    np.arange(np.float16(0.1), np.float16(0.5), np.float16(0.1)).dtype,
    np.arange(np.float32(0.1), np.float32(0.5), np.float32(0.1)).dtype,
    np.arange(np.float64(0.1), np.float64(0.5), np.float64(0.1)).dtype,
    np.arange(np.float128(0.1), np.float128(0.5), np.float128(0.1)).dtype
    ])

print([
    np.r_[np.float16(0.1):np.float16(0.5):np.float16(0.1)].dtype,
    np.r_[np.float32(0.1):np.float32(0.5):np.float32(0.1)].dtype,
    np.r_[np.float64(0.1):np.float64(0.5):np.float64(0.1)].dtype,
    np.r_[np.float128(0.1):np.float128(0.5):np.float128(0.1)].dtype
    ])

print([
    np.mgrid[np.float16(0.1):np.float16(0.5):np.float16(0.1)].dtype,
    np.mgrid[np.float32(0.1):np.float32(0.5):np.float32(0.1)].dtype,
    np.mgrid[np.float64(0.1):np.float64(0.5):np.float64(0.1)].dtype,
    np.mgrid[np.float128(0.1):np.float128(0.5):np.float128(0.1)].dtype
    ])

print([ # note comma
    np.mgrid[np.float16(0.1):np.float16(0.5):np.float16(0.1),].dtype,
    np.mgrid[np.float32(0.1):np.float32(0.5):np.float32(0.1),].dtype,
    np.mgrid[np.float64(0.1):np.float64(0.5):np.float64(0.1),].dtype,
    np.mgrid[np.float128(0.1):np.float128(0.5):np.float128(0.1),].dtype
    ])

Output:

# input   float16      float32           float64          float128
[dtype('float64'), dtype('float64'), dtype('float64'), dtype('float128')] # linspace, logspace, geomspace
[dtype('float64'), dtype('float64'), dtype('float64'), dtype('float128')] # arange
[dtype('float64'), dtype('float64'), dtype('float64'), dtype('float128')] # r_
[dtype('float64'), dtype('float64'), dtype('float64'), dtype('float128')] # mgrid
[dtype('float64'), dtype('float64'), dtype('float64'), dtype('float64')]  # mgrid w/ comma

The mgrid w/ comma only runs with the recent bug fix (#16815), in which I thought all floats in these similar methods were cast to float64, so it is inconsistent with the rest. @eric-wieser mentioned making mgrid infer lower precision types could be follow up work, which prompted this issue. I'm happy to fix that, but I'm not sure what the correct output should be.

Personally I think the code for all these functions should be updated to match what the docs imply unless there is a reason the output needs the precision of at least a float64. Other option would be to just fix mgrid w/ comma to match the rest and maybe further clarify the docs.

Numpy/Python version information:

1.20.0.dev0+f457a1a 3.7.7 (default, Mar 10 2020, 15:43:33)
[Clang 11.0.0 (clang-1100.0.33.17)]

Start by reviewing the dtype behavior of linspace, geomspace, logspace, arange, r_, and mgrid, including the comma form of mgrid described in the report. Determine the intended lower-precision behavior, then make the affected entry points consistent and update the corresponding documentation and regression coverage.