Let $U$ be a bounded domain of $\mathbb{R}^d$,  and write $m$ for the Lebesgue measure on $U$. For $k=1,2$, we denote by $H^k(U)$ the set of all locally $m$-integrable functions $u\colon U \to \mathbb{R}$
such that for any multi-index $\alpha$ with $|\alpha|\le k$, the weak derivative $D^\alpha u$ exists and belongs to $L^2(U,m)$.

Define the Neumann Laplacian $(L,\text{Dom}(L))$ on $U$ by
\begin{align*}
\text{Dom}(L)&=\left\{u \in L^2(U,m) :  H^1(U) \ni v \mapsto \int_{U}\nabla u\cdot \nabla v\,dm \text{ is continuous on $L^2(U,m)$}\right\} \\
-\int_{U}v Lu\,dm&=\int_{U} \nabla u\cdot \nabla v\,dm,\quad u \in \text{Dom}(L),\,v \in H^1(U).
\end{align*}

I know that if $U$ is a bounded $C^2$ domain, $\text{Dom}(L) \subset H^2(U)$ (see Section~10.6.2 in [[1]] ). Even if $U$ is a bounded Lipschitz domain, does this inclusion hold? I don't think this is correct in general. 

What is the most general condition for $U$ such that this inclusion holds?


  [1]: https://link.springer.com/book/10.1007/978-94-007-4753-1