What does $e^{Ax}$ mean? Consider the differential equation $$y'=Ay+b$$ , where $A$ is a $n\times n$-matrix,
 y and b are vectors of functions $(y_1(x),...y_n(x))^T$ and $(b_1(x),...,b_n(x))^T$
Suppose, we have found the general solution of the homogenous equation $$y'=Ay$$
The standard approach is to continue with the "variation of the constants". But in 
 a script, I found that the matrix $e^{Ax}$ has something to do with the special solution
 of the inhomogenous equation.
What is the meaning of $e^{Ax}$ ?
For example, let $$A=\pmatrix{2&-1\\-1&2}$$ , $$b=\pmatrix{-e^{3x}\\e^{3x}}$$
What is $e^{Ax}$ and a special solution of $$y'=Ay+b$$ ?
 A: It is the matrix exponential:
$$
e^{Ax}=\sum_{n=0}^\infty\frac{x^n}{n!}\,A^n.
$$
It is a fundamental solution of the homogeneous equation $y'=A\,y$, and is the unique one equal to the identity matrix when $x=0$. The unique solution such that $y(0)=y_0$ is $y(x)=e^{Ax}y_0$.
The method of variation of constants looks for a particular solution of the complete equation of the form $e^{Ax}C(c)$, where $C(x)=(C_1(x),\dots,C_n(x))$.
A: in your specific problem, it is easy to find the eigenvalues and eigenvectors of the coefficient matrix $\pmatrix{2&-1\\-1&2}.$ the eigenvalues are $3$ and $1$ the corresponding orthonormal eigenvectors are first and second columns of the matrix $V = \pmatrix{1/\sqrt 2& 1/\sqrt 2\\-1/\sqrt 2& 1/\sqrt 2}.$ we have $$A =  V \pmatrix{3&0\\0&1}V^T, \, e^{Ax} =   V \pmatrix{e^{3x}&0\\0&e^x}V^T = \frac{1}{2}\pmatrix{e^{x} + e^{3x}& e^x - e^{3x}\\ e^x - e^{3x}&e^{x} + e^{3x}}$$
a particular solution of the nonhomogeneous problem $\frac{dy}{dx} = Ay + b$ is 
$$y_p = \int_0^x e^{A(x-\xi)}e^{3\xi}\pmatrix{-1\\1}\, d\xi=\int_0^x e^{3\xi} \pmatrix{0\\e^{3(x-\xi)}}\, d\xi = e^{3x}\pmatrix{0\\x}.$$
the general solution is $$y = e^{Ax}y(0) + y_p. $$
A: There are several ways to solve this problem. One way consists in applying the Laplace transform to the differential equation so that:
\begin{equation}
L\{y'(x)\}=A L\{y(x)\}+L\{b(x)\}
\end{equation}
By writing $L\{y(x)\}=Y(s)$ and $L\{b(x)\}=B(s)$ the equation above reduces to:
\begin{equation}
s Y(s)-y(0)= A Y(s) + B(s)
\end{equation}
Then
\begin{equation}
Y(s) (I s -A) = I y(0) + B(s)
\end{equation}
And
\begin{equation}
Y(s)  = y(0) (s I-A)^{-1} + B(s) (s I-A)^{-1}
\end{equation}
And finally we go back in the x domain by applying the inverse Laplace transform
\begin{equation}
y(x)  = y(0) L^{-1}\{(s I-A)^{-1}\} + b(x) \star L^{-1}\{ (s I-A)^{-1}\}
\end{equation}
So that when comparing the above relation to the solution of the differential equation
\begin{equation}
y(x) = y(0) e^{A x}+\int_0^x e^{A(x-\psi)} b(\psi) d\psi
\end{equation}
we can conclude that
\begin{equation}
e^{A x} = L^{-1}\{(s I-A)^{-1}\}
\end{equation}
