Monge-Kantorovich duality with a $\{0,1\}$ cost function Consider the usual Monge-Kantorovich transportation problem where $X$ and $Y$ are Polish spaces, $\mu$ and $\nu$ are probability measures on $X$ and $Y$, and $c:X\times Y \to \mathbb{R}^+ \cup \{+\infty \}$ is a lower semi-continuous cost function.  The Kantorovich duality theorem states that the transportation cost between $\mu$ and $\nu$ is equal to the supremum of $$\int_X \varphi~ d\mu +\int_Y \psi~ d\nu $$ over all $L_1$ functions $\varphi(x)$ and $\psi(y)$  such that  $\varphi(x)+\psi(y)\leq c(x,y)$ for $d \mu$-almost all $x\in X$ and $d \nu$-almost all $y\in Y$.
My question is:  if $c(x,y)\in \{0,1\}$ for all $x$ and $y$, and for each $x$ there exists $y_1,y_2$ such that $c(x,y_1)=0$ and $c(x,y_2)=1$ (and similarly for each $y$ there exist $x_1,x_2$ such that the equivalent condition holds) does it follow that there exists a solution (or "almost exists" a solution) where $\varphi(x)$ and $\psi(x)$ only take values in the set $\{-1,0,1\}$?  Finite dimensional experiments with linear programs suggest that the answer is "yes" but I cannot tell if they extend to the general setting.
UPDATE:  I added additional conditions about $c$, which guarantees that we can find solutions such that $0\leq\varphi\leq1$ and $-1\leq \psi \leq 0$.  This is because we can shift $\varphi$ and $\psi$ such that $\sup_x \varphi(x)=1$, which guarantees that $\psi(y)\geq-1$.  Furthermore, it must always be true that $\psi(y)\leq 0$ because otherwise there would exist $x,y$ such that $\varphi(x) + \psi(y) > 1$.  This in turn guarantees that we can always assume that $\varphi(x)\geq 0$ as well.
 A: The strong Kantorovich duality, i.e. the existence of dual solutions, holds whenever $c$ is lower semicontinuous and real-valued, and the optimal coupling has finite cost, see Theorem 5.10(ii) in Villani "Optimal Transport old and new". In contrast, for the weak duality ("$\inf = \sup$") it suffices to assume $c$ is lower semicontinuous.
To prove 5.10(ii) one constructs an explicit $c$-concave function $\varphi$ that is finite-valued on the support of $\mu$, see formula (5.17). Beware of signs! Villani is looking at $c$-convex functions and you at $c$-concave functions.
To answer the "value part" of the question: the values of $\psi$ and $\varphi$ are only important on the support of $\mu$ and resp. $\nu$. Everywhere else you may choose them to be $-\infty$ (loosing $c$-concavity). Now some cases: If the optimal coupling $\pi$ has cost $1$ then for all $(x,y),(x',y') \in \Gamma := \operatorname{supp}\pi$ it holds $c(x,y')=1$ as otherwise $c(x',y)+c(x,y') < 2 = c(x,y)+c(x',y')$ violating $c$-cyclic monotonicty. Thus $\varphi \equiv 0 \equiv \psi-1$ is a dual solution to the problem. A similar case happens if the cost is $0$ on $\operatorname{supp}\mu \times \operatorname{supp}\nu$. For the last case, pick $(x_0,y_0) \in \Gamma$ with $c(x_0,y_0)=0$ and do the standard construction to get an integer-valued $c$-concave $\varphi$ with $\varphi(x_0)=0$. Observe if $(x_1,y_1),(x_2,y_2)\in \Gamma$ then   
\begin{align*}
\varphi(x_{1}) + \psi(x_1) & = c(x_{1},x_{1})\\
\varphi(x_{2}) + \psi(x_2) & = c(x_{2},x_{2}).
\end{align*}
Subtracting this from $\varphi+\psi\le c$ applied to the couples $(x_1.y_2)$ and $(x_2,y_1)$ gives
\begin{align*}
\varphi(x_{2})-\varphi(x_{1}) & \le c(x_{2},y_{1})-c(x_{1},x_{1})\\
\varphi(x_{1})-\varphi(x_{2}) & \le c(x_{1},y_{2})-c(x_{2},y_{2}).
\end{align*}
The right hand side has values in $\{-1,0,1\}$. 
Since $\varphi(x_0)=0$ the function $\varphi$ has values in the same set. Since $1-(-1)=2$, it must have values either in $\{0,1\}$ or in $\{-1,0\}$. A similar argument as above then shows that $\psi$ needs have values in $\{-1,0\}$ or $\{0,1\}$. 
Note if $c$ is not real-valued then there may not be any dual solutions between certain measures $\mu$ and $\nu$, see the preprint of S. Suhr and me on Lorentzian cost function, more explicitly Section 3.1.
