  The Bose--Hubbard dimer model is a celebrated fundamental quantum mechanical
model that accounts for the dynamics of bosons at two interacting sites. It has
been realized experimentally by two coupled, driven and lossy photonic crystal
nanocavities, which are optical devices that operate with only a few hundred
photons due to their extremely small size. Our work focuses on characterizing
the different dynamics that arise in the semiclassical approximation of such
driven-dissipative photonic Bose--Hubbard dimers. Mathematically, this system
is a four-dimensional autonomous vector field that describes two specific
coupled oscillators, where both the amplitude and the phase are important. We
perform a bifurcation analysis of this system to identify regions of different
behavior as the pump power $f$ and the detuning $\delta$ of the driving signal
are varied, for the case of fixed positive coupling. The bifurcation diagram in
the $(f,\delta)$-plane is organized by two points of codimension-two
bifurcations -- a $Z_2$-equivariant homoclinic flip bifurcation and a Bykov
T-point -- and provides a roadmap for the observable dynamics, including
different types of chaotic behavior. To illustrate the overall structure and
different accumulation processes of bifurcation curves and associated regions,
our bifurcation analysis is complemented by the computation of kneading
invariants and of maximum Lyapunov exponents in the $(f,\delta)$-plane. The
bifurcation diagram displays a menagerie of dynamical behavior and offers
insights into the theory of global bifurcations in a four-dimensional phase
space, including novel bifurcation phenomena such as degenerate singular
heteroclinic cycles.
