Yes there is. Here is how you do it. First find an orthogonal change in variables $$ x_j=\sum_{jk} s_{jk}y_k $$ $(s_{jk})$ orthogonal matrix, so that in the new coordinates we have $$ \sum_{i,j}a_{ij} x_ix_j = \sum_j \mu_j^2 y_j^2, $$ where $\mu_j^2$ are the eigenvalues of the symmetric matrix $(a_{ij})$. Note that $\newcommand{\ii}{\boldsymbol{i}}$ $\newcommand{\pa}{\partial}$ $$\frac{\pa}{\pa x_k}=\sum_j\frac{\pa y_j}{\pa x_k}\frac{\pa }{\pa y_j} =\sum_j t_{jk}\frac{\pa}{\pa y_j}, $$ where $(t_{jk})$ is the inverse of the orthogonal matrix $(s_{jk})$ so that $t_{jk}=s_{kj}$. Then for some real numbers $r_j$ $$ -\sum_j\frac{\pa^2}{\pa x_j^2}+2\ii\sum_j\lambda_j\frac{\pa }{\pa x_j}+\sum_{i,j}a_{ij}x_ix_j $$ $$ =-\sum_j\frac{\pa^2}{\pa y_j^2} +2\ii\sum_j r_j\frac{\pa}{\pa y_j} +\sum_j \mu_j^2 y_j^2 $$ $$=\underbrace{\sum_j \left(\ii\frac{\pa}{\pa y_j}+r_j\right)^2 +\sum_j\mu_j^2y_j^2-\sum_j r_j^2}_{=: L}. $$ Next set $$ R^2 :=\sum_j r_j^2, \;\; w(t,y) := R^2t +\sum_j \ii r_j y_j, $$ $$ L_0 :=\sum_j\left(-\frac{\pa^2}{\pa y_j^2} +\mu_j^2y_j^2\right), $$ and observe that $$\pa_t +L =e^{w}(\pa_t+L_0) e^{-w}. $$ Thus, if $K$ is a fundamental solution of $\pa_t+L_0$, then $$(\pa_t +L) (e^{w} K) = e^{w} (\pa_t+L_0)K=e^{w} \delta_0=\delta_0 $$ so that $e^{w}K$ is a fundamental solution of $\pa_t+L$. You can find $K$ using Mehler's formula.