Problem

subplots() defaults to auto-squeezing the returned Axes array (returning a single array when one is requested, returning a 1D array when a single row or single column is requested), likely for practicality, but this makes dynamic layouts pretty annoying: if you write axs = fig.subplots(m, n) where m and n are somehow dynamically computed (based on some dataset characteristics) and then later plan to use axs[i, j], you basically have an IndexError waiting to happen unless you can be certain that m and n are never 1. Of course one can remember to always write subplots(..., squeeze=False), but that's a bit verbose.

Proposed Solution

My preferred solution would probably have been to make the default m and n not 1, but None, so that None means "1 in that direction, and squeeze that direction" whereas 1 just means 1 (without squeezing). This would have fixed most common cases, e.g. if you write subplots() you get a single axes, subplots(3) or subplots(nrows=3) or subplots(ncols=3) you get a 1D array, subplots(m, n) (where m and n are dynamically computed variables) you always get a 2D array. I think in practice the main back-incompatibility of this proposal that occurs in real life is subplots(1, 3) (single row, likely not so uncommon), which would typically have to be written subplots(None, 3) or subplots(ncols=3) or fig, (axs,) = subplots(3) (for the additional unpacking).

Assuming that this backcompat breakage is not acceptable (it likely is not), one alternative would be to change the signature of subplots() to also support taking a single 2-tuple as argument (rather than nrows and ncols separately, and make subplots((m, n)) never squeeze the input. (Likely we wouldn't bother supporting m or n being None, in that case.) Then it would be OK to retrain my muscle memory to just always add an extra pair of parentheses when calling subplots.

Additional context and prior art