Let $u\in H^1(\Omega)$ be the weak solution of the Laplacian equation in the unit cube $\Omega \subset \mathbb{R}^n$, $n \leq 3$, with inhomogeneous Neumann conditions on $\Gamma_N$ and homogeneous Dirichlet conditions on $\Gamma_D$. That is, $u$ solves $$- \Delta u = 0 \ \text{ in } \Omega, \qquad \frac{\partial u}{\partial n}= g \ \text{ on } \Gamma_N, \qquad u = 0 \ \text{ on } \Gamma_D.$$
Fix an arbitrary open ball $B \subset \Omega$. Given $c \in \mathbb{R}^n$, is it possible to choose Neumann data $g$ (in a suitable space) such that $$\nabla u = c \quad \text{ a.e. on } B?$$
If this is not the case, then how can I characterize attainable gradient fields?
Idea: Green's function might be a suitable tool to address this question.