  We give a new and rigorous duality relation between two central notions of
weak solutions of nonlinear PDEs: entropy and viscosity solutions. It takes the
form of the nonlinear dual inequality: \begin{equation}\int |S_t u_0-S_t v_0|
\varphi_0 \mathrm{d}x\leq \int |u_0-v_0| G_t \varphi_0 \mathrm{d}x, \quad
\forall \varphi_0 \geq 0, \forall u_0, \forall v_0, \qquad(\star)\end{equation}
where $S_t$ is the entropy solution semigroup of the anisotropic degenerate
parabolic equation \begin{equation*} \partial_t u+\mathrm{div} F(u) =
\mathrm{div} (A(u) D u),\end{equation*} and where we look for the smallest
semigroup $G_t$ satisfying ($\star$). This amounts to finding an optimal
weighted $L^1$ contraction estimate for $S_t$. Our main result is that $G_t$ is
the viscosity solution semigroup of the Hamilton-Jacobi-Bellman
equation\begin{equation*} \partial_t \varphi = \mathrm{sup}_\xi \{F'(\xi) \cdot
D \varphi+\mathrm{tr}(A(\xi) D^2\varphi)\}.\end{equation*} Since weighted $L^1$
contraction results are mainly used for possibly nonintegrable $L^\infty$
solutions $u$, the natural spaces behind this duality are $L^\infty$ for $S_t$
and $L^1$ for $G_t$. We therefore develop a corresponding $L^1$ theory for
viscosity solutions $\varphi$. But $L^1$ itself is too large for
well-posedness, and we rigorously identify the weakest $L^1$ type Banach
setting where we can have it -- a subspace of $L^1$ called
$L^\infty_{\mathrm{int}}$. A consequence of our results is a new domain of
dependence like estimate for second order anisotropic degenerate parabolic
PDEs. It is given in terms of a stochastic target problem and extends in a
natural way recent results for first order hyperbolic PDEs by [N. Pogodaev, J.
Differ. Equ., 2018].
