LTI systems in convolution representation LTI systems in state space representation are systems of the form:
\begin{eqnarray}
\dot{\mathbf{x}}(t)=\mathbf{Ax}(t)+\mathbf{Bu}(t)
\end{eqnarray}
\begin{eqnarray}
\mathbf{y}(t)&=&\mathbf{Cx}(t)+\mathbf{Du}(t)
\end{eqnarray}
where $\mathbf{x}$ is the state vector, $\mathbf{u}$ is the input and $\mathbf{y}$ is the output. These systems satisfy the superposition principle only if the initial condition $\mathbf{x}=\mathbf{0}$ holds. Let us refer to the class of these systems with the aforementioned initial condition as $\Sigma_{LTI}^0$.
Now a different class of systems, namely $\Sigma_{conv}$ are the ones that are represented by a convolution integral, i.e. the input-output behaviour of the system is described by:
\begin{equation}
\mathbf{y}(t)=\int_{-\infty}^\infty h(\tau) \mathbf{u}(t-\tau)d\tau
\end{equation}
My question is whether it is true that every system in $\Sigma_{LTI}^0$ admits a representation in $\Sigma_{conv}$? If yes, how?
Note: I understand that my question boils down to  finding a function $h$ such that (assume for simplicity that $\mathbf{C}=\mathbf{I}$ and $\mathbf{D}=\mathbf{0}$):
\begin{equation}
\int_0^t e^{\mathbf{A}(t-\tau)}\mathbf{Bu}(\tau)d\tau= \int_{-\infty}^\infty h(\tau) \mathbf{u}(t-\tau)d\tau
\end{equation}
Update: Let us assume that there exists a $h\in\mathcal{L}^2(\Re;\Re^n)$ such that an LTI system is equivalent to a convolution system and assume for simplicity that the input is one-dimensional. Then, these have identical impulse responses, hence for $u(t)=\delta(t)$ we have:
\begin{equation}
\int_0^t e^{\mathbf{A}(t-\tau)}\mathbf{B}\delta(\tau)d\tau= \int_{-\infty}^\infty h(\tau) \delta(t-\tau)d\tau
\end{equation}
Therefore:
\begin{equation}
\int_0^t e^{\mathbf{A}(t-\tau)}\mathbf{B}\delta(\tau)d\tau= h(t) \implies
h(t)=e^{\mathbf{A}t}\mathbf{B}
\end{equation}
But then:
\begin{equation}
\int_0^t e^{\mathbf{A}(t-\tau)}\mathbf{Bu}(\tau)d\tau= \int_{-\infty}^\infty e^{\mathbf{A}(t-\tau)}\mathbf{Bu}(\tau)d\tau
\end{equation}
Doesn't look good! So, if there is no mistake in the above procedure, there is no $h$ that fulfills my requirements. Then, what is the connection between the aforementioned classes of systems.
 A: Initial conditions can be represented as impulse responses applied to the system at time $t=0$. In other words, you can add a new column $x(0)$ to the $B$ matrix and extend your input function to $\pmatrix{u(t) \\ \delta(t)}$. In a way, you treat the initial condition as a perturbation to the zero initial condition system since 
$$
x(t) = Ce^{At}x_0 + \int_0^t Ce^{\mathbf{A}(t-\tau)}\mathbf{Bu}(\tau)d\tau
$$
Note that, $\delta(t)$ function only makes sense under the integral sign.
Therefore, every LTI system in $\Sigma^0_{LTI}$ is a member of $\Sigma^{x_0}_{LTI}$ with the special case $x_0=0$.
A: I think your problem is that your forgot to add that your solution for $h(t)=e^{\mathbf{A}t}\mathbf{B}$ is valid for $t>0$ 
And that makes sense: assuming for one moment ${\bf D}=0$, we have a system in which the output is given by a first-order linear differential equation, which indeed corresponds to an exponential $h(t)$ (causal, i.e., $h(t) =0$ for $t<0$). Adding ${\bf D}\ne0$ you're just adding a ${\bf D} \; \delta(t)$ term, so the conclusion would be that $\Sigma_{LTI}^0$ is a (small) subset of $\Sigma_{conv}$, namely, that of the LTI which $h(t)$ is given by the sum of an (casual) exponential and a Dirac delta.  
This reasoning should be extended to multidimensional systems, though. Now, assume ${\bf u}$ and ${\bf x}$ are one-dimension, but we are allowed to increase the dimension of ${\bf x}$ arbritrarily (which is equivalent to extend arbitrarily the order of the ODE). Then, it's another story.
BTW, if you are familiar with discrete LTI systems, the analogy would be $\Sigma_{LTI}^0 \leftrightarrow AR(1)$ (autoregresive processes of order 1, i.e., one pole), which indeed corresponds to a subset of all LTI casual systems. 
A: The transfer function of the given system is
$$
\mathbf{Y}(s)=[\mathbf{D}+\mathbf{C}(s\mathbf{I}-\mathbf{A})^{-1}\mathbf{B}]\mathbf{U}(s)
$$
and it establishes an input-output relationship written compactly as
$$
\mathbf{Y}(s)=\mathbf{G}(s)\mathbf{U}(s)
$$
Using the convolution property of the Laplace transform, we have:
$$
\mathbf{y}(t)=\mathbf{g}(t)*\mathbf{u}(t)
$$
where $g(t)=(\mathscr{L}^{-1}\mathbf{G})(t)$. Therefore, every system in state-space form can be written in an equivalent input-output representation using the convolution operation.
