Let $\rho\in C(\bar{D})$ be a continuous function on the compact set $\bar{D}$ and let us define $$\mathbf{E}(\mathbf{x}_0):=k\int_D\frac{\rho(\mathbf{x})}{\|\mathbf{x}_0-\mathbf{x}\|^3}(\mathbf{x}_0-\mathbf{x})d\mathbf{x}$$where I have used the short notation $d\mathbf{x}$ for $dxdydz$, with $\mathbf{x}=(x,y,z)$, and where $k$ is a constant, Coulomb's constant if we intend $\rho$ to be an electric charge density and $\mathbf{E}$ the electrostatic field. Clearly $$\int_D\frac{\rho(\mathbf{x})}{\|\mathbf{x}_0-\mathbf{x}\|^3}(\mathbf{x}_0-\mathbf{x})d\mathbf{x}=-\int_{D-\mathbf{x_0}}\frac{\rho(\mathbf{x}+\mathbf{x}_0)}{\|\mathbf{x}\|^3}\mathbf{x}d\mathbf{x}$$therefore I think that imposing conditions upon $\rho$ would allow us to have a finite $\mathbf{E}$ variously subject to desired conditions of regularity.
Although this problem arises in a physics context, I would like to find a mathematical proof of how to guarantee the usual conditions of regularity assumed in physics for $\mathbf{E}$. For example if $\rho\in C^k(A)$ where $A$ is open and contains $\bar{D}$, can it be guaranteed that $\mathbf{E}$ is of class $C^k$? If it can, how is it proved?
The fact that the domain of integration $D-\mathbf{x}_0$ depends upon the variable(s) $\mathbf{x}_0$ does not allow me to use standard results of differentiation under the integral sign like the fact that if $V\subset\mathbb{R}^3$ is compact and $f:V\times[a,b]\to\mathbb{R}$ has a continuous partial derivative $\frac{\partial f}{\partial t}\in C(V\times[a,b])$ then for all $t\in[a,b]$ $$\frac{d}{dt}\int_V f(x,y,z, t)dxdydz=\int_V\frac{\partial}{\partial t} f(x,y,z, t)dxdydz.$$
I heartily thank you for any answer!