Changing the relative phase of two photons in an Hong-Ou-Mandel interference produces anti-bunching? What happens to the output of a beam splitter when you change the relative phase between two photons that enter from the two input ports?
In Hong-Ou-Mandel interference for a beamsplitter of the form, where I represent my outputs as $b^\dagger_1$ and $b^\dagger_2$:
$$ \begin{equation*}
 \left(\begin{array}{cc}
  \hat{b}^\dagger_1\\
\hat{b}^{\dagger}_2   \\
 \end{array}\right) = \frac{1}{\sqrt{2}}
\left(\begin{array}{cc}
 1  &   -1 \\ 
 1 &  1 \\
 \end{array}\right)
 \left(\begin{array}{cc}
 a^\dagger_1  \\
 a^\dagger_2   \\
 \end{array}\right) 
\end{equation*} $$
which implies the inputs have the relation:
$$ \begin{equation*}
 \left(\begin{array}{cc}
  \hat{a}^\dagger_1\\
\hat{a}^{\dagger}_2   \\
 \end{array}\right) = \frac{1}{\sqrt{2}}
\left(\begin{array}{cc}
 1  &   1 \\ 
 -1 &  1 \\
 \end{array}\right)
 \left(\begin{array}{cc}
 b^\dagger_1  \\
 b^\dagger_2   \\
 \end{array}\right) 
\end{equation*} $$
with an input of $|1, 1\rangle = a^\dagger_1 a^\dagger_2 |0, 0\rangle = \frac{1}{\sqrt{2}}(b^\dagger_1+b^\dagger_2)\frac{1}{\sqrt{2}}(-b^\dagger_1+b^\dagger_2)= \frac{1}{2}(-b^\dagger_1 b^\dagger_1-b^\dagger_2b^\dagger_1 +b^\dagger_1 b^\dagger_2+b^\dagger_2 b^\dagger_2) = \frac{1}{2}(-b^\dagger_1 b^\dagger_1+b^\dagger_2 b^\dagger_2)$
This math, to me, suggests that the resultant "two-photon interference" is invariant to the relative phase between the two fields. That is, if I add a phase $e^{i \phi}$ to one of my $a^\dagger$ modes, it just gets carried through the whole process as a global phase, without producing interference:
$|\tilde{1}, 1\rangle = \left(a^\dagger_1 e^{i \theta}\right) a^\dagger_2 |0, 0\rangle = e^{i \theta}(-b^\dagger_1 b^\dagger_1+b^\dagger_2 b^\dagger_2)$
This phase doesn't change the fact that the photons $|1, 1\rangle$ states destructively interfere. I thought this fact is aligned with the general intuition that "photons don't have well-defined phases'' because generally pure Fock states will often lose any phase given to them unless a relative phase is created (for instance putting a Fock state in a Mach-Zehnder interferometer).
But this conclusion appears to be in contradiction with this paper, which says that adding a relative phase to the photon pair ends up changing the interference, allowing to flip between bunching and antibunching depending on the phase.
In this paper they say that you can think of the result as a sort of post-selected Mach-Zehnder interferometer. If the first photon is found in detector 1, it means the second photon acts like it is in a Mach-Zehnder interferometer, and consequently can be routed into either detector changing the relative phase between paths.  To quote:

So what exactly is wrong about this previous treatment?
 A: I don't know if this will help, as I don't fully understand the way the paper is modelling the situation, but here's how I would describe it.
The point is that we want to study what happens to the many-body interference when the photons stop being indistinguishable. In this case, they stop being (completely) indistinguishable due to differences in their wavefunctions.
In other words, the two-photon state that is injected into the interferometer has the form
$a^\dagger_{A,\psi}a^\dagger_{B,\phi},$ where $|\psi\rangle$ and $|\phi\rangle$ are the photons' wavefunctions, and $A,B$ denote the two inputs of the interferometer (and we omit the vacuum state on which these operators act for notational brevity). If the two photons are indistinguishable, then $\psi=\phi$, and the state can be written more simply as $a^\dagger_A a^\dagger_B$. Evolving this through the beamsplitter gives you the usual HOM effect etc.
However, what if the photons' wavefunctions are not identical, or the photons are distinguishable by any other mean (e.g. by their time of arrival to the interferometer)? We can generally describe this type of situation by writing
$$|\phi\rangle = \alpha|\psi\rangle + \beta|\psi_\perp\rangle,$$
where $|\alpha|^2+|\beta|^2=1$ and $|\psi_\perp\rangle$ is some state orthogonal to $|\psi\rangle$. The input state is then
$$a^\dagger_{A,\psi}(\alpha \,a^\dagger_{B,\psi}+\beta \,a^\dagger_{B,\psi_\perp})
= \alpha\, a^\dagger_{A,\psi} a^\dagger_{B,\psi} + \beta \,a^\dagger_{A,\psi}a^\dagger_{B,\psi_\perp}.$$
Evolving through the interferometer, the first term gives you the usual HOM, while the second term behaves like two distinguishable photons. By tuning the overlap between $|\phi\rangle$ and $|\psi\rangle$, i.e. by tuning $\alpha$ and $\beta$, we can see the transition between distinguishability and indistinguishability.
A: Let’s try to work backwards.
$$\begin{aligned}
|1, 1\rangle &= a^\dagger_1a^\dagger_2 |0, 0\rangle \\
&= 
\frac{b^\dagger_1+b^\dagger_2}{\sqrt{2}}~\frac{-b^\dagger_1+b^\dagger_2}{\sqrt{2}} |0, 0\rangle\\
&= 
\frac{1}{2}(-b^\dagger_1 b^\dagger_1-b^\dagger_2b^\dagger_1 +b^\dagger_1 b^\dagger_2+b^\dagger_2 b^\dagger_2) |0, 0\rangle
\end{aligned}$$
The only way we’ll get an output different from bunching is if the cross terms here have a phase difference so that they don’t cancel out. The only way this can happen is if they pick up a different phase between the corresponding $b_i^\dagger$s that originate from different $a_j^\dagger$s. That is to say:
$$\frac{b^\dagger_1+b^\dagger_2}{\sqrt{2}}~\frac{-b^\dagger_1+b^\dagger_2}{\sqrt{2}} \to  \frac{e^{i\phi}b^\dagger_1+b^\dagger_2}{\sqrt{2}}~\frac{-b^\dagger_1+b^\dagger_2}{\sqrt{2}} $$
So now the question is how do we physically realise such a thing? Remember that we get bunching only when the photons are indistinguishable. So this hints us towards photons being distinguishable. In fact, the paper you’re referring to  does it in a clever way. They distinguish it in time. For one of the photons they apply a phase shift of $\phi$ for half the wavepacket by applying an appropriate step voltage in time (at the peak of the packet). For more details look at figure 2 in the paper.

So for the first half there will be no coincidence count. But there will be in the second half depending on the voltage applied. For more details look at figure 3 of the paper.
A: I suspect the point is that one can introduce an additional phase shifter so that your scattering matrix is now equivalent to a matrix of the form
$$ \begin{equation*}
 \left(\begin{array}{cc}
  \hat{a}^\dagger_1\\
\hat{a}^{\dagger}_2   \\
 \end{array}\right) = \frac{1}{\sqrt{2}}
\left(\begin{array}{cc}
 1  &   e^{i\phi} \\ 
 -1 &  1 \\
 \end{array}\right)
 \left(\begin{array}{cc}
 b^\dagger_1  \\
 b^\dagger_2   \\
 \end{array}\right) 
\end{equation*} $$
While I don't know the experimental details to produce such a shift I know its possible to obtain this as the proposed matrix is $\in U(2)$.  Moreover, this is still a 50/50 device since the transmittance and reflectance are the mod-squared of the entries, and still clearly both equal to $1/2$.
Anyways this then produces
\begin{align}
a_1^\dagger a_2^\dagger \to \frac{1}{2}
\left(-b_1^\dagger b_1^\dagger + e^{i\phi} b_2^\dagger b_2^\dagger - e^{i\phi}b_1^\dagger b_2^\dagger + b_2^\dagger b_1^\dagger\right)
\end{align}
Detecting one photon in each port is then possible using $\hat\Pi=\vert 11\rangle\langle 11\vert$ with probability
$$
P_{11}(\phi)= \frac{1}{2}\left(1-\cos(\phi)\right)=\sin^2(\phi/2) \, . \tag{A}
$$
In particular for $\phi=0$ we recover the HOM result for indistinguishable photons.
I suppose this is equivalent to the proposed 1-photon states $\vert\tilde\Psi_{\pm} \rangle$ in the sense that these 1-photon states yield the same prob. as (A).
A: The detection of the first photon projects the input state into the state of a photon which is in the superposition state of both input modes:
$\Psi_{det} = 1\sqrt{2}(\pm|1\rangle_2 +|1\rangle_1)$
where the relative phase is given by the detector in which it is detected. If you now take this state and propagate it through the beamsplitter, you obtain exactly what you describe, namely that the second photon will always be detected in the same output port. However, if you change the phase of the input state (and this is what is done with the EOM), e.g.
$\Psi_{det} = 1\sqrt{2}(\pm|1\rangle_2 + \exp(i \Delta_{\phi})|1\rangle_1)$
you can also change the output port. This is basically the same as in a single photon interference experiment where you can change in which output port is detected by changing the phase of the arms in the interferometer.
Edit: In the paper that you linked (the one with the theory) they interfere two photons with a slightly different frequency. By then changing the time interval between the detection events, one can also change the output port. It's section three in the paper (if Im not mistaken). That's based on the same effect.
