Solution to system of linear ODEs backwards in time Suppose we have a system of linear ODEs represented in state-space form as follows:
$$ \dot{x}(t) = Ax(t) + Bu(t), \; x(t_0)=x_0. $$
The well-known solution, valid for $t\in[t_0,\infty)$, is
$$ x(t) = e^{A(t-t_0)}x_0 + \int_{t_0}^t e^{A(t-\sigma)}Bu(\sigma)d\sigma. $$
What is the general form of the solution when $x(t_f)=x_f$ is specified instead, and we want to know the solution for $t\in(-\infty,t_f]$?

If we suppose that $u$ is defined in an interval $(a,b)$( can be $\mathbb R$ ) and $t_0\in (a,b)$, then the solution of
$$ \dot{x}(t) = Ax(t) + Bu(t), \; x(t_0)=x_0. $$
is
$$ x(t) = e^{A(t-t_0)}x_0 + \int_{t_0}^t e^{A(t_0-s)}Bu(s)ds. $$
FOR ALL $t\in (a,b)$, not only $t\geq t_0$.
To see this 
write the equation like this 
$$x'-Ax=Bu(t)$$
Now we apply $e^{(t_0-t)A}$ in the equation
$$e^{(t_0-t)A}(x'-Ax)=e^{(t_0-t)A}Bu(t)\implies \dfrac{d}{dt}e^{(t_0-t)A}x=e^{(t_0-t)A}Bu(t)$$
Therefore 
$$\int_{t_0}^{t}\dfrac{d}{ds}e^{(t_0-s)A}xds=\int_{t_0}^te^{(t_0-s)A}Bu(s)ds$$
and we get 
$$ x(t) = e^{A(t-t_0)}x_0 + \int_{t_0}^t e^{A(t_0-s)}Bu(s)ds. $$
