Is the relationship $\vec{D} = \epsilon \vec{E}$ always valid? In electromagnetism class we had a problem proposed which led me to some thoughts about the relationship between $\vec{D}$ and $\vec{E}$. I know it certainly $\textbf{is}$ a definition so it must be true always but let me explain:
In this excercise we have a spherical capacitor (inner radius $r_1 < r_2$ outer radius) with a dielectric permittivity $\varepsilon = \varepsilon_0 \left( k_1 + k_2 \cos^2 \theta \right)$, which is obviously dependent on the spherical coordinate $\theta$. A conductor is an equipotential surface so the electric potential $V$ can't depend on $\theta$. If we assume the inner terminal is positive charged with $Q$, we can compute: as $\vec{\nabla}\cdot \vec{D} = \rho$, we have $\vec{D} = \frac{Q}{4\pi r^2}$ and $Q = \int \int_S \vec{D} \cdot \hat{n} \mathrm{d}S'$ = $\int \int_S \varepsilon \left(\theta\right) E(r) \mathrm{d}S' = 4 \pi \varepsilon_0 E(r) r^2 \left( k_2 + k_2/3 \right)$, so $\vec{E} = \frac{Q}{4\pi \varepsilon_0 r^2} \frac{1}{\left( k_2 + k_2/3 \right)} \neq \varepsilon\left(\theta \right)^{-1} \vec{D}$. I know this makes sense, as because of the simmetry of this problem, the displacement field cannot depend on $\theta$ and, because of the fact that $V$ doesn't depend on $\theta$, $\vec{E}$ has only a radial component and also independent on $\theta$. 
What is happening? Why is $\vec{D} \neq \varepsilon(\theta) \vec{E}$?
 A: Whilst I think you can assume that the electric field is only radial and with a magnitude that depends only on radius (because of the equipotentials and the requirement that the E-field is normal to the conducting surfaces), I don't think you can assume the same for the D-field.
If we assume that $\vec{E} = E_r(r) \hat{r}$ and a linear dielctric, then the integral form of Gauss's law applied to a spherical surface of radius $R$ passing through the dielectric, tells us
$$Q = 
\oint \vec{D} \cdot d\vec{A} = \int^{2\pi}_0 \int^{\pi}_0 \epsilon(\theta) E_r(R) R^2 \sin \theta\ d\theta\ d\phi$$
In this case $\epsilon(\theta) = \epsilon_0 (k_1 + k_2\cos\theta)$. So
$$Q = 2\pi \epsilon_0 E_r(R)R^2 \int^{\pi}_0 k_1 \sin \theta + \frac{1}{2}k_2\sin(2\theta)\ d\theta  = 4\pi k_1 \epsilon_0 E_r(R)R^2   $$
Thus 
$$\vec{E} = \frac{Q}{4\pi k_1 \epsilon_0 r^2}\ \hat{r}$$
Assuming a linear, isotropic medium, then 
$$\vec{D} = \epsilon \vec{E} = \frac{(k_1 + k_2 \cos \theta)Q}{4\pi k_1 r^2}\ \hat{r}$$
If you integrate this over a spherical surface, you will of course find that it comes to $Q$, but the fact that $\vec{D}$ depends on $\theta$ is telling you that the charge is not uniformly distributed over the surfaces.
If we zoom in on a piece of the inner sphere (so close that it looks flat) and draw a Guassian pillbox straddling the inner conductor surface, then we can say that, since there is no D-field inside the conductor, then if the charge per unit area is $\sigma$,
$$\vec{D}(r_1) \cdot d\vec{A} = \sigma \ dA$$ 
i.e.
$$ \sigma = D(r_1) = \epsilon E = \frac{(k_1 + k_2 \cos \theta)Q}{4\pi k_1 r_1^2}$$
Thus I think the root cause of your difficulties is assuming that the D-field has no $\theta$ dependence and that the charge is uniformly distributed on the conductors.
As to what circumstances you can assume $\vec{D} = \epsilon \vec{E}$, though that is a red herring in the case of this question, it is that the medium must be linear and isotropic. Note the comment by @jamison , you can still use the linear relationship at each point even in an inhomogeneous medium. 
A: $\mathbf D=\epsilon \mathbf E$ is only valid for when there material is a linear dielectric that is isotropic and homogeneous. In your case this is not valid with your $\theta$ dependence, so this is why you do not see this relationship.
Remember, $\mathbf D$ is defined as 
$$\mathbf D=\epsilon_0\mathbf E+\mathbf P$$
When you have the properties described above, you then have $\mathbf P=\epsilon_0\chi\mathbf E$ which gets you to what you were proposing
$$\mathbf D=\epsilon_0\mathbf E+\epsilon_0\chi\mathbf E=\epsilon\mathbf E$$
A: Since $\vec D$ is a response of a medium to the applied electric field it can differ along each direction i.e the response is not isotropic. Thus in general $\vec D$ is given by the following:
$$\vec D =\epsilon_0\overleftrightarrow \chi \vec E$$
where $\overleftrightarrow \chi$ is the susceptibility tensor. This is still under the linear response regime. For non-linear regime you may have complicated dependence of $\vec D$ on $\vec E$.
