ldiv for vector-of-vectors RHS

[stevengj]

Close #904 to allow the right-hand-side for ldiv calls to be a vector of vectors: i.e. solving $Ax=b$ where the elements of $b$ and $x$ are themselves vectors (e.g. arrays, or elements of some arbitrary vector field).

Mostly, this was a straightforward matter of changing oneunit to zero (vector fields should support zero but need not have one), with two caveats:

• This works as long as the eltype itself supports zero, e.g. it works for Vector{<:SVector} but not Vector{<:Vector} because zero(Vector{T}) fails (it doesn't know the size). This seems like a reasonable limitation to me.
• In one place, for ldiv(F::Factorization, B::AbstractVecOrMat), for some reason the code first converts the factorization F itself to a type TFB based on promotion with B. I don't actually understand why this conversion is necessary — it dates back to Make adjoint and solve work for most factorizations julia#40899 by @andreasnoack. Andreas, any comments? If we keep that conversion, however, then it needs to know the "scalar" type of eltype(B). We don't currently provide a generic function for this. In this PR, I make my best guess via a new function _scalartype that calls eltype recursively until it hits a Number type, or until it hits a type with EltypeUnknown. (This should at least be backwards-compatible?) I would rather just remove the FF = Factorization{TFB}(F) conversion entirely, if possible, since it makes this code less generic.

(A workaround to the latter issue is to call ldiv!, since the in-place routines don't convert the factorization type.)

To do:

• Decide what to do about FF = Factorization{TFB}(F) — remove this conversion, keep the _scalartype workaround, or do something else?
• Tests — I tested this locally using b = [@SVector rand(2) for _=1:2]; A = rand(2,2); x = A \ b; @test A*x ≈ b, but to avoid a StaticArrays dependency we may need to add a dummy SVector-like type to test/testhelpers.jl