Trotter expansions in ode solver?

Question

Trotter expansions say: $$ e^{A+B} = \lim_{P\to\infty} \big(e^{A/2P} e^{B/P} e^{A/2P} \big)^P. $$ With $P = 2$, it becomes (with high accuracy) $$ e^{A/4} e^{B/2} e^{A/2} e^{B/2} e^{A/4}. $$

Let's say I want to solve $\frac{dy}{dt} = Cy$, where $C = A + B$. $C$ can be time dependent. How can I solve this equation by splitting $C$ into $A$ and $B$, before using any ode solving scheme. Will it yield a gain in performance by doing so, like in matrix exponential?

Background information: https://www.quora.com/What-is-the-significance-of-the-Trotter-product-formula-in-physics

Answer 1

What you're proposing is widely referred to as Strang splitting. There is a huge literature on similar methods. There are also many questions on this site on the same topic and a tag for them.

Your question seems to suggest that in your case the ODE is linear (but possibly time dependent), and you don't assume any special structure of $C$. In this case there is (generically) nothing to be gained by using a splitting method. Splitting methods are useful when your right-hand-side can be decomposed into two additive parts, each of which is much easier to integrate than the combination.

Answer 2

Your equation $\frac{dy}{dt} = C(t)y$ has the exact solution

$$ y(t) = \mathcal T exp\left(\int_0^t C(s) ds\right) y(0) $$

where $\mathcal T$ is the time-ordering operator which orders operators at later time to the left. Approximate the integral by a discrete sum ($t=N\Delta t$),

$$y(t) = \mathcal T exp\left(\sum_{i=1}^N H(t_i) \Delta t\right) y(0)$$

and split the terms (*),

$$y(t) = \mathcal T \prod_{i=1}^N exp\left(H(t_i) \Delta t\right) y(0).$$

Now you can trivially apply the time-ordering operator, leading to

$$y(t) = \prod_{i=1}^N exp\left(H(t_i) \Delta t\right) y(0).$$

The latter equation says that one can repeatedly apply matrix exponentials to obtain the solution. This is a fact on which almost any numerical integrator relies.

Now note the step marked by (*). Here, the exponential of a sum of several operators is decomposed into a product of these exponentials, that is, the elementary splitting formula $exp(H_1+H_2) = exp(H_1) exp(H_2)$ is applied several times (the formula is, of course, only an approximation when the two operators do not commute).

How does this help you in solving your equation? In no way. But it's good to know. For the practical part of your question, please read the first answer in your linked question on quora, it basically says it all.