I am interested in finding references regarding estimates of the form $$ \| D^2 u\|_{L^2(\Omega)} \leq C(\|f\|_{L^2(\Omega)}+\|g\|_{L^2(\Gamma_N)} )$$ where $\|D^2 u\|_{L^2(\Omega)}^2 = \sum\limits_{i,j \in \{1,2\}} \int_{\Omega}(\partial_{x_i x_j}u)^2$ and $$ \begin{cases} -\Delta u = f & \text{ in }\Omega \\ u = 0 & \text{ on } \Gamma_D \\ \frac{\partial u}{\partial n} = g & \text{ on } \Gamma_N \end{cases} $$ and $\Omega\subset \Bbb{R}^2$ is polygonal (bounded, convex) such that the solution $u$ belongs to $H^2(\Omega)$ (in particular $g$ is regular enough). The boundaries $\Gamma_D, \Gamma_N$ are both non-void and make a partition for $\partial\Omega$.

I was able to find in the book of Grisvard (Elliptic Problems in Nonsmooth Domains) a similar result when $g\equiv 0$. However, I could not find this precise estimate for the non-homogeneous Neumann boundary conditions.

Can you point me to a reference that contains proofs of such estimates? In particular, I am interested in finding precisely the constant $C$ (if possible).