Bilinear embedding for perturbed divergence-form operator with complex coefficients on irregular domains

Let $\Omega\subseteq\mathbb{R}^{d}$ be open, $A$ a complex uniformly strictly accretive $d\times d$ matrix-valued function on $\Omega$ with $L^{\infty}$ coefficients, $b$ and $c$ two $d$-dimensional vector-valued functions on $\Omega$ with $L^{\infty}$ coefficients and $V$ a locally integrable nonegative function on $\Omega$. Consider the operator ${\mathscr L}^{A,b,c,V}=-{\rm div}\,(A\nabla) + \left\langle \nabla , b \right\rangle - {\rm div}\,(c \, \cdot) + V$ with mixed boundary conditions on $\Omega$. We extend the bilinear inequality that Carbonaro and Dragi\v{c}evi\'c proved in the special cases when $b=c = 0$. As a consequence, we obtain that the solution to the parabolic problem $u^{\prime}(t)+{\mathscr L}^{A,b,c,V}u(t)=f(t)$, $u(0)=0$, has maximal regularity in $L^{p}(\Omega)$, for all $p>1$ such that $A$ satisfies the $p$-ellipticity condition that Carbonaro and Dragi\v{c}evi\'c introduced in arXiv:1611.00653 and $b,c,V$ satisfy another condition that we introduce in this paper. Roughly speaking, $V$ has to be ``big'' with respect to $b$ and $c$. We do not impose any conditions on $\Omega$, in particular, we do not assume any regularity of $\partial\Omega$, nor the existence of a Sobolev embedding.