[proofplan]
We remove the homogeneous dynamics by multiplying by $e^{-tA}$. The transformed inhomogeneous term is an indefinite [Lebesgue integral](/page/Lebesgue%20Integral) of an $L^1$ function, hence is absolutely continuous and differentiates to its integrand almost everywhere. Multiplying back by $e^{tA}$ gives the stated formula and verifies the differential equation. Uniqueness follows by applying the same integrating-factor transformation to the difference of two trajectories.
[/proofplan]
[step:Define the candidate trajectory by an integrating factor]
Choose a Borel measurable representative of the equivalence class $u \in L^1([0,T];\mathbb{R}^m)$, and denote this representative again by
\begin{align*}
u: [0,T] \to \mathbb{R}^m.
\end{align*}
For a [linear map](/page/Linear%20Map) $L: \mathbb{R}^m \to \mathbb{R}^n$, let
\begin{align*}
\|L\|_{\mathrm{op}}=\sup\{|Lz|: z \in \mathbb{R}^m, |z|\le 1\}
\end{align*}
denote its Euclidean operator norm. Define the measurable map
\begin{align*}
f: [0,T] \to \mathbb{R}^n,\qquad f(\tau)=e^{-\tau A}Bu(\tau).
\end{align*}
Since $\tau \mapsto e^{-\tau A}B$ is continuous on the compact interval $[0,T]$, there is a finite constant
\begin{align*}
M=\sup_{\tau \in [0,T]}\|e^{-\tau A}B\|_{\mathrm{op}}<\infty.
\end{align*}
Thus $|f(\tau)| \le M|u(\tau)|$ for $\mathcal{L}^1$-a.e. $\tau \in [0,T]$, so $f \in L^1([0,T];\mathbb{R}^n)$.
Define the indefinite integral
\begin{align*}
F: [0,T] \to \mathbb{R}^n,\qquad F(t)=\int_0^t f(\tau)\,d\mathcal{L}^1(\tau).
\end{align*}
By the absolute [continuity theorem](/theorems/1145) for indefinite Lebesgue integrals applied to the function $f \in L^1([0,T];\mathbb{R}^n)$, $F$ is absolutely continuous and
\begin{align*}
\dot{F}(t)=f(t)=e^{-tA}Bu(t)
\end{align*}
for $\mathcal{L}^1$-a.e. $t \in [0,T]$.
Now define
\begin{align*}
x: [0,T] \to \mathbb{R}^n,\qquad x(t)=e^{tA}(x_0+F(t)).
\end{align*}
The map $t \mapsto e^{tA}$ is continuously differentiable, and $F$ is absolutely continuous, so $x$ is absolutely continuous.
[guided]
The point of the integrating factor is to convert the forced equation into an equation whose right-hand side is directly integrable. Since $u \in L^1([0,T];\mathbb{R}^m)$ is an equivalence class of [measurable functions](/page/Measurable%20Functions) modulo equality almost everywhere, we first choose a Borel measurable representative and denote it again by
\begin{align*}
u: [0,T] \to \mathbb{R}^m.
\end{align*}
Changing this representative on a set of $\mathcal{L}^1$-measure zero does not change any of the Lebesgue integrals below. For a linear map $L: \mathbb{R}^m \to \mathbb{R}^n$, define its Euclidean operator norm by
\begin{align*}
\|L\|_{\mathrm{op}}=\sup\{|Lz|: z \in \mathbb{R}^m, |z|\le 1\}.
\end{align*}
We now define
\begin{align*}
f: [0,T] \to \mathbb{R}^n,\qquad f(\tau)=e^{-\tau A}Bu(\tau).
\end{align*}
This is the forcing term after removing the homogeneous flow. We must check that this expression is integrable before using it inside an indefinite integral. The map $\tau \mapsto e^{-\tau A}B$ is continuous because the matrix exponential is smooth in $\tau$, and $[0,T]$ is compact. Therefore
\begin{align*}
M=\sup_{\tau \in [0,T]}\|e^{-\tau A}B\|_{\mathrm{op}}
\end{align*}
is finite. For $\mathcal{L}^1$-a.e. $\tau \in [0,T]$,
\begin{align*}
|f(\tau)|=|e^{-\tau A}Bu(\tau)| \le \|e^{-\tau A}B\|_{\mathrm{op}}|u(\tau)| \le M|u(\tau)|.
\end{align*}
Since $u \in L^1([0,T];\mathbb{R}^m)$, this domination proves $f \in L^1([0,T];\mathbb{R}^n)$.
We may therefore define
\begin{align*}
F: [0,T] \to \mathbb{R}^n,\qquad F(t)=\int_0^t f(\tau)\,d\mathcal{L}^1(\tau).
\end{align*}
The absolute continuity theorem for indefinite Lebesgue integrals applies because $f \in L^1([0,T];\mathbb{R}^n)$, and it gives that $F$ is absolutely continuous with derivative
\begin{align*}
\dot{F}(t)=f(t)=e^{-tA}Bu(t)
\end{align*}
for $\mathcal{L}^1$-a.e. $t \in [0,T]$.
Finally define
\begin{align*}
x: [0,T] \to \mathbb{R}^n,\qquad x(t)=e^{tA}(x_0+F(t)).
\end{align*}
This is exactly the result of solving for the transformed variable and then multiplying back by the homogeneous flow. Because $t \mapsto e^{tA}$ is continuously differentiable and $F$ is absolutely continuous, their product is absolutely continuous.
[/guided]
[/step]
[step:Rewrite the candidate in the stated convolution form]
For every $t \in [0,T]$, the definition of $F$ gives
\begin{align*}
x(t)=e^{tA}x_0+e^{tA}\int_0^t e^{-\tau A}Bu(\tau)\,d\mathcal{L}^1(\tau).
\end{align*}
Since $e^{tA}$ is independent of the integration variable $\tau$, linearity of the Lebesgue integral gives
\begin{align*}
x(t)=e^{tA}x_0+\int_0^t e^{tA}e^{-\tau A}Bu(\tau)\,d\mathcal{L}^1(\tau).
\end{align*}
Because $A$ commutes with itself, the matrix exponentials satisfy $e^{tA}e^{-\tau A}=e^{(t-\tau)A}$. Hence
\begin{align*}
x(t)=e^{tA}x_0+\int_0^t e^{(t-\tau)A}Bu(\tau)\,d\mathcal{L}^1(\tau).
\end{align*}
[/step]
[step:Verify the initial condition and differential equation]
At $t=0$, the integral over $[0,0]$ is zero, so
\begin{align*}
x(0)=e^{0A}x_0=x_0.
\end{align*}
For $\mathcal{L}^1$-a.e. $t \in [0,T]$ at which $\dot{F}(t)=e^{-tA}Bu(t)$, the product rule for absolutely continuous functions gives
\begin{align*}
\dot{x}(t)=Ae^{tA}(x_0+F(t))+e^{tA}\dot{F}(t).
\end{align*}
Substituting the derivative of $F$ yields
\begin{align*}
\dot{x}(t)=Ae^{tA}(x_0+F(t))+e^{tA}e^{-tA}Bu(t).
\end{align*}
Let $I_n \in \mathbb{R}^{n \times n}$ denote the identity matrix. Since $x(t)=e^{tA}(x_0+F(t))$ and $e^{tA}e^{-tA}=I_n$, this becomes
\begin{align*}
\dot{x}(t)=Ax(t)+Bu(t).
\end{align*}
Thus the candidate is an absolutely continuous solution of the initial value problem.
[/step]
[step:Prove uniqueness by differentiating the transformed difference]
Let $y: [0,T] \to \mathbb{R}^n$ be another absolutely continuous map satisfying $y(0)=x_0$ and
\begin{align*}
\dot{y}(t)=Ay(t)+Bu(t)
\end{align*}
for $\mathcal{L}^1$-a.e. $t \in [0,T]$. Define
\begin{align*}
w: [0,T] \to \mathbb{R}^n,\qquad w(t)=e^{-tA}(y(t)-x(t)).
\end{align*}
The map $w$ is absolutely continuous. For $\mathcal{L}^1$-a.e. $t \in [0,T]$ where both $x$ and $y$ satisfy the differential equation, the product rule gives
\begin{align*}
\dot{w}(t)=-Ae^{-tA}(y(t)-x(t))+e^{-tA}(\dot{y}(t)-\dot{x}(t)).
\end{align*}
Since
\begin{align*}
\dot{y}(t)-\dot{x}(t)=A(y(t)-x(t)),
\end{align*}
we obtain
\begin{align*}
\dot{w}(t)=-Ae^{-tA}(y(t)-x(t))+e^{-tA}A(y(t)-x(t)).
\end{align*}
The matrices $A$ and $e^{-tA}$ commute, so the two terms cancel and
\begin{align*}
\dot{w}(t)=0
\end{align*}
for $\mathcal{L}^1$-a.e. $t \in [0,T]$. An absolutely [continuous function](/page/Continuous%20Function) with derivative zero almost everywhere is constant, hence $w(t)=w(0)$ for every $t \in [0,T]$. Since
\begin{align*}
w(0)=e^{0A}(y(0)-x(0))=0,
\end{align*}
we have $w(t)=0$ for every $t \in [0,T]$. Multiplying by $e^{tA}$ gives $y(t)=x(t)$ for every $t \in [0,T]$. Therefore the trajectory given by the formula is unique.
[/step]