Time-stepping Krylov: special case for expv_timestep

[MSeeker1340]

As suggested by @ChrisRackauckas, added a special case for phiv_timestep which calculates just the first term $exp(tA)*b$.

We can already do this using phiv_timestep(!) alone, but that requires using a n-by-1 matrix B, which is a bit unintuitive. expv_timestep(!) is basically a wrapper for phiv_timestep(!) that allows a vector b as input. (Internally it is reshaped to the n-by-1 matrix form).

I was initially a bit worried about doing extra work for the special case. After more carefully examining the code for updating the W matrix and u, I found that all the loops are no-ops for p == 0, so we're not losing performance here.

[ChrisRackauckas]

I was initially a bit worried about doing extra work for the special case. After more carefully examining the code for updating the W matrix and u, I found that all the loops are no-ops for p == 0, so we're not losing performance here.

Great! Yes, most things should be no-ops. The only thing to be careful about is the interior expm! call. The Sidje algorithm in the case of just expv has the inner part simplify to just using exp, right? I assume that right now you're doing expm! on a 1x1 in expv_timestep? If that's the case, that's something to specialize on.

[MSeeker1340]

The Sidje algorithm in the case of just expv has the inner part simplify to just using exp, right?

I see the misunderstanding here. The Sidje algorithm computes expm! on a (m + p) x (m + p) matrix where m is the Krylov subspace size and p is the highest order phi. Even for p=0 (that is, expv_timestep) we still need to compute expm! on the original m x m Heisenberg matrix H. So there's nothing to optimize/specialize here.

[ChrisRackauckas]

oh perfect. Then yeah, this is good.