Consider the multivariate Ornstein--Uhlenbeck process (https://en.wikipedia.org/wiki/Ornstein%E2%80%93Uhlenbeck_process):
$$ dx_t = - \beta x_t dt + \sigma dW_t$$
where $\beta$ and $\sigma$ are constant $d \times d$ matrices. I would like to solve the OU process by first solving for the mean and then for the deviation of the mean. According to wiki, given a deterministic $x_0$, the mean is simply $\mu(x_t) = \exp(-\beta t) x_0$. The difference between mean and the SDE should then be $x_t - \mu(x_t) = \exp(-\beta) \sigma \int_0^t \exp(t-t') dW_{t'}$. The last missing piece is how to compute $\int_0^t \exp(t-t') dW_{t'}$. Is this analytically solvable?
It seems that in the univariate case $\int_0^t f(t') dW_{t'} = \mathcal{N}(0, \int_0^t f(t')^2 dt') $(https://quant.stackexchange.com/questions/57066/integration-over-function-of-wiener-process). However, I have a multivariate Brownian motion. Does the result simply extend to a multivariate normal distribution with zero mean and covariance matrix with all diagonal entries $\int_0^t f(t')^2 dt'$ and off-diagonal entries being 0?