Integrate gaussian times sqare root of x times polynomial of order 2
As in your previous question, we have the following results, in terms of Bessel functions,
depending on whether $\mu$ is positive or negative. But first, let us define the following three
quantities : $\color{green}A=\dfrac\mu{2\sigma}~,~$ $\color{green}P=\mu~\Big[2~\big(a-\mu\big)^2+3~\sigma^2\Big],~$ and $\color{green}Q=2\mu~\big(a-\mu\big)^2+\big(5\mu-4a\big)~\sigma^2.~$
Then, for $\color{blue}{\mu<0}$ we have $F~=~\dfrac{e^{-A^2}}{8\sigma^2}\cdot\sqrt{-\dfrac\mu2}~\cdot~\bigg[Q\cdot K_{\tfrac14}\Big(A^2\Big)-P\cdot K_{\tfrac34}\Big(A^2\Big)\bigg],~$ and
for $\color{blue}{\mu>0}~$ we get $$F~=~\dfrac\pi2\cdot\dfrac{e^{-A^2}}{8\sigma^2}\cdot\sqrt{\mu}~\cdot~\bigg\{P\cdot\bigg[I_{\tfrac34}\Big(A^2\Big)+I_{-\tfrac34}\Big(A^2\Big)\bigg]+Q\cdot\bigg[I_{\tfrac14}\Big(A^2\Big)+I_{-\tfrac14}\Big(A^2\Big)\bigg]\bigg\},~$$ whereas for $\color{blue}{\mu=0}$ we obtain $F~=~\dfrac{\big(2a^2+3\sigma^2\big)~\Gamma\bigg(\dfrac34\bigg)-a~|\sigma|~\sqrt2\cdot\Gamma\bigg(\dfrac14\bigg)}{4~\sqrt[\large^4]2~\sqrt{|\sigma|}}~.$
There is, surprisingly, an easy way to do this using sums. Firstly, all I actually have to do is evaluate integrals of the form $$ I(s,m) = 2\int_0^{\infty} y^{s-1} \exp{(-(y-m)^2)} \, dy: $$ your integral can be expressed in terms of a sum of integrals of this form, by scaling $y$ and substituting $m=\mu/(\sqrt{2}\sigma)$ and so on. So, first I expand the exponential and extract the $m^2$ term, leaving $$ e^{m^2}I(s,m) = 2\int_0^{\infty} y^{s-1} e^{-y^2} e^{-2my} \, dy. $$ There is an obvious next step: expand the right-hand side as a power series in $m$. Interchanging the order of integration, we have $$ e^{m^2}I(s,m) = 2\int_0^{\infty} y^{s-1} e^{-y^2} \left( \sum_{k=0}^{\infty} \frac{(-1)^k}{k!} (2m)^k y^k \right) \, dy \\ = \sum_{k=0}^{\infty} \frac{(-1)^k}{k!} (2m)^k \left( 2\int_0^{\infty} y^{k+s-1} e^{-y^2} \, dy \right) $$ The integral can be done by substituting $u=y^2$, $du/u=2dy/y$, which gives $$ 2\int_0^{\infty} y^{k+s-1} e^{-y^2} \, dy = \int_0^{\infty} u^{(k+s)/2-1} e^{-u} \, du = \Gamma\left(\tfrac{1}{2}(k+s)\right)$$
It remains to do the sum, $$ e^{m^2}I(s,m) = 2\int_0^{\infty} y^{s-1} e^{-y^2} \left( \sum_{k=0}^{\infty} \frac{(-1)^k}{k!} (2m)^k y^k \right) \, dy \\ = \sum_{k=0}^{\infty} \frac{(-1)^k}{k!}\Gamma\left(\tfrac{1}{2}(k+s)\right) (2m)^k $$ It is straightforward to turn this into a couple of confluent hypergeometric functions by separating the even and odd terms and playing with the $\Gamma$s using the duplication formula (I'll have a go if asked), to reach the final result $$ I(s,m) = e^{-m^2} \left( \Gamma\left(\tfrac{1}{2}s\right) {}_1 F_1\left(\tfrac{1}{2}s;\tfrac{1}{2}; m^2 \right) -2m \Gamma\left(\tfrac{1}{2}(s+1)\right) {}_1 F_1\left(\tfrac{1}{2}(s+1);\tfrac{3}{2}; m^2 \right) \right), $$ where $$ {}_1 F_1 (a;b;z) = \frac{\Gamma(b)}{\Gamma(a)}\sum_{k=0}^{\infty} \frac{\Gamma(a+k)}{\Gamma(b+k)} \frac{z^k}{k!}. $$