Integral of $e^{-at}$ Erf[t] and $e^{-at} Erf[t]^2$

Can anyone help me with any suggestions for the following triple integral?

$$\frac{1}{\sqrt{2\pi \sigma^2}}\int_{-\infty}^t dt_1 e^{-(\gamma_1+\gamma_2)t_1} \int_{-\infty}^{t_1} dt_2 e^{\frac{\gamma_1+\gamma_2}{2}t_2} e^{-\frac{t_2^2}{4\sigma^2}} \int_{-\infty}^{t_2} dt_3 e^{\frac{\gamma_1+\gamma_2}{2}t_3} e^{-\frac{t_3^2}{4\sigma^2}}$$ It can be put into the equivalent form, $$\frac{1}{2\sqrt{2\pi \sigma^2}}\int_{-\infty}^t dt_1 \left[ \left[ e^{-\frac{\gamma_1+\gamma_2}{2}t_1} \Theta(t_1) \right] * e^{-\frac{t_1^2}{4\sigma^2}} \right]^2$$ where $$\Theta(t_1)$$ is the Heaviside unit step function and * denotes convolution.

This second integrand can be written in terms of he error function, Erf[],

$$e^{\frac{(\gamma_1+\gamma_2)^2 \sigma^2}{2}} \frac{1}{2\sqrt{2\pi \sigma^2}}\int_{-\infty}^t dt_1 e^{-(\gamma_1+\gamma_2)t_1} \left[ \frac{1}{2} + \frac{1}{2}Erf\left[\frac{t_1-(\gamma_1+\gamma_2) \sigma_1^2}{2\sigma} \right] \right]^2$$ The definition of Erf[] being $$Erf[\frac{t-t_0}{\sigma}] = \int_0^{t-t_0} dt_1 \frac{1}{\sqrt{2\pi\sigma^2}} e^{-\frac{t_1^2}{2\sigma^2}}$$

Mathematica fails with the computation, even numerically, in the last part.

As asked in title $$I=\int e^{-a t}\, \text{erf}(t)\,dt$$ is simple (just one integration by parts and completing the square) $$I=\frac{e^{\frac{a^2}{4}}\, \text{erf}\left(t+\frac{a}{2}\right)-e^{-a t}\, \text{erf}(t)}{a}$$
But the second one $$J=\int e^{-a t}\, \big[\text{erf}(t)\big]^2\,dt$$ seems to be much more problematic.