Exotic Topology in Ordinary Soft Matter

Hydrodynamics is shown to induce non-Hermitian topological phenomena in ordinary, passive soft matter. This is demonstrated for the first time by subjecting a 2D elastic lattice to a low-Reynolds viscous flow. The interplay of hydrodynamics and elasticity splits Dirac cones into bulk Fermi arcs, pairing exceptional points with opposite half-integer topological charges. The bulk Fermi arc is a generic hallmark of the system exhibited in all lattice and flow symmetries. Analytic model and simulations explain how the emergent singularities shape the spectral bands and give rise to a web of van Hove singularity lines in the density of states. The present findings suggest that non-Hermitian physics can be explored in a broad class of ordinary soft matter, living and artificial alike, opening avenues for topology-based technology in this regime.


INTRODUCTION
The conservation of energy in isolated Hermitian systems is a basic tenet of physics, but in practice, most systems are open, exchanging energy and information with the external world. This inherent non-Hermiticity is traditionally seen as an inevitable imperfection of realistic systems, yet recent studies revealed that it gives rise to distinctive phenomena unmatched in Hermitian physics [1]-most notably skewed spectral bands prone to symmetry breaking when exceptional points emerge [2][3][4][5][6]. And this discovery kicked off intensive efforts to engineer non-Hermitian systems in diverse classical and quantum settings, ranging from photonics [7][8][9][10], phononics [11][12][13], and optomechanics [14] to electronics [15] and atomic lattices [16]. The potential prowess of non-Hermitian technology has been demonstrated in developing new meta-materials and devices [17][18][19][20].
But an open question remains: besides sophisticated engineered systems, can one observe and utilize these so-called exotic topological phenomena in more common and natural settings?-After all, we are immersed in a dissipative, non-Hermitian world, and living systems are immanently open to exchange with the surrounding environment at all scales [21]. The present study offers a clear positive answer: simple model and simulations provide demonstrate for the first time non-Hermitian topological hallmarks in standard, passive elastic networks subject to ordinary viscous flow [22,23]. Such settings are omnipresent in the overdamped low-Reynolds regime, typical to cells, macromolecules and simple soft matter systems, suggesting that non-Hermitian topology is much more common than previously realized.
Among the exotic phenomena observed in non-Hermitian materials, bulk Fermi arcs [7] hold a special place. In contrast to the ingrained intuition that frequency levels are closed curves, each Fermi arc is an open isofrequency curve ending at two exceptional points. * tsvitlusty@gmail.com These endpoints are defects with opposite half-integer topological charges. It is important to note that bulk Fermi arcs are topological hallmarks of non-Hermiticity in the bulk spectrum of the lattice, unlike the more familiar surface Fermi arcs induced by Weyl points in 3D Hermitian systems [24][25][26][27]. Bulk arcs have been so far elusive and were empirically observed in one photonic crystal [7]. Recent theoretical studies predict bulk arcs in the spectra of heavy fermions [28] and spin liquids [29]. Strikingly, the present work finds that the elusive bulk arcs are a generic phenomenon of driven viscoelastic matter [30][31][32][33][34] whose emergence does not require fine tuning and stems directly from the symmetry of the interactions.
An important step into the dissipative regime are recent theoretical studies [35,36] demonstrating non-Hermitian "odd" elastodynamics [37] in a network of active units specifically designed to exert circular forces. Likewise, non-Hermitian viscosity waves were shown to emerge in a system with "odd viscosity" [38]. Common to these "odd" systems are nonreciprocal elasticity or viscosity coefficients, originating from active modules that are hard to realize in microscopic settings. In sharp contrast, the present system requires no active parts, but merely standard passive soft matter in a laminar background flow, where non-reciprocity is a direct outcome of the hydrodynamic interaction [22,31]. Thus, linking hydrodynamics and topology in soft matter can expand the realm of non-Hermitian physics to yet further fields, beyond the already rich plethora of applications developed in recent years.
In the following Results section, we derive the dynamic equations of the hydro-elastic lattice and map its spectrum to a standard non-Hermitian form in 2D momentum space. Then, we demonstrate how the interplay of Hermitian and skew-Hermitian interactions generates bulk Fermi arcs, which are shown to be a generic hallmark of the flow-induced dynamics. Next, we discuss the topological features of these singularities and their links to the density of states. Finally, in the Discussion, we suggest possible experimental realizations and examine potential implications and future directions.  1. (A) A triangular lattice made of spheres connected by elastic struts is immersed in a thin layer of viscous fluid between two walls in the z direction (perpendicular to the page). The lattice is driven at a velocity u relative to the fluid. (B) A pair of moving particles induce dipolar flow fields (grey streamlines), thereby exerting on each other equal hydrodynamic forces, f hyd ij = f hyd ji (blue arrows) in the direction of the induced flow. The elastic forces along the elastic strut are opposite, f elas ij = − f elas ji (red arrows), thus conserving linear momentum.

RESULTS
The motion of the hydro-elastic lattice. To see how non-Hermitian topology arises in ordinary elastic matter at low-Reynolds, consider the following model system (Fig. 1A). A two-dimensional triangular lattice made of spherical particles of size -joined by thin elastic struts of length a and spring constant κ-is submerged in a viscous fluid of viscosity η and is moving at a velocity u relative to the fluid. The resulting viscous drag on each particle is γu, where γ is the friction coefficient. Such relative motion can be obtained in the lab by holding the lattice in a flow or by or by driving the lattice with various forces (see Discussion).
Perturbing the surrounding fluid, the driven particles induce long-range interactions throughout the lattice [22,[30][31][32] (Fig. 1B). When the flow is limited to a thin sheet of viscous fluid between solid floor and ceiling, these hydrodynamics forces are dipolar [39,40], and the force exerted by particle j on particle i is where R i are the particle positions and the distance vec- The magnitude of the hydrodynamic dipole in (1) scales as ∼u 2 . The Hookean elastic forces are proportional to the change of the length, ∆R ij = R ij −R ij , whereR ij the equilibrium length of the springs (a in the lattice), Owing to the dipolar symmetry of (1), the hydrodynamic forces a pair of particles exert on each other are equal, f hyd ij = f hyd ji (because θ ji = π + θ ij ). Thus, the dipolar forces break Newton's third law of momentum conservation (Fig. 1B). This is because viscous flow is an inherently open system, an effective representation of energy and momentum transfer from hydrodynamic degrees-of-freedom to microscopic ones. To steadily move the lattice, the momentum leakage ('loss') needs to be constantly compensated by the driving force ('gain'). In contrast, the elastic forces (2) are opposite and conserve linear momentum, f elas ij = −f elas ji . As shown below, the interplay of conservative and non-conservative forces give rise to skewed non-Hermitian topology.
Two timescales govern the dynamics of the lattice: the hydrodynamic timescale, τ hyd = a 3 /(u 2 ), and the elastic relaxation time, τ elas = γ/κ. Their ratio is the hydroelastic number, which controls the system's behavior: when 1, it is dominated by hydrodynamics, and when 1 by elasticity. Hereafter, we measure times in units of the hydrodynamic timescale τ hyd and frequencies in τ −1 hyd . The dynamical equations combine Stokes flow with Hookean elasticity, a linear regime where the correspondence between experiment and theory is well-established. In this overdamped regime, the friction force on each particle is counterbalanced by the driving force, F = Fx, and by the hydrodynamic and elastic interactions in the lattice, whereṘ i is the i th particle's velocity (see App. A 1). In equation (4), the long-range hydrodynamic forces (1) are summed over all particles, while elastic interactions (2) are summed only among neighbors connected by struts. Owing to the lattice parity symmetry, the sums of interactions vanish at steady-state, when the lattice traverses uniformly at a velocity u = F/γ.
Dynamics in momentum space. The coaction of elastic and hydrodynamic forces excites collective motion in the lattice. Expanding the dynamics (4) in small x + Ω 2 y . The spectrum exhibits six Dirac points (green, one denoted as D) on the boundary of the Brillouin zone (black hexagon). Middle: At each Dirac point, the bands merge, forming a graphene-like double-cone, "diabolo" shape. Right: the 3D shape of the frequency bands in the first Brillouin zone (grey hexagon), showing the double-cones, which are halved by the zone's boundary. (B) At = 0, the symmetry is broken when H includes a skew-Hermitian component, H elas . Left: The real part of the bands, Re(ω+) = − Re(ω−), drawn for = π, exhibits six bulk Fermi arcs (white lines) in the first Brillouin zone (the whole spectrum is shown in Fig. 7). Each arc splits from a Dirac point (green) and joins two Exceptional Points (ExPs) with topological charges ± 1 2 (orange and light blue, E+ and E− are two ExPs that split from D). Middle: the real parts of the bands merge along the Fermi arc, forming a double-wedge shape. Right: 3D shape of the the frequency bands in the first Brillouin zone (grey hexagon), showing the double-wedges (halved by the zone's boundary).
deviations r j of the particles from their steady-state positions in the moving latticeR j , we find that the collective modes are plane waves r j = e k exp i(k ·R j − ωt) (App. A 2-A 3). The 2D polarization of the wave, e k , is an eigenstate of a Schrödinger-like equation with an eigenfrequency ω, where the operator H is a momentum-space representation of the forces in the lattice. The "Hamiltonian" H is a 2×2-matrix, which sums the hydrodynamic and elastic interactions with a relative weight . In the basis of left-and right-circular polarizations, we obtain and where σ x , σ y and 1 are Pauli's and the unity matrices. The frequencies Ω x , Ω y , ω x , ω y and ω 1 , are the Fourier sums of the interactions-all real by the parity symmetry of the lattice (see App. A 3). The hydroelastic operator H (6) is analogous to the Hamiltonian of spin-1 2 particles in a complex magnetic field with damping [41], and this spinor-like nature shows in the spectrum, as discussed below. Equation (5) with its non-Hermitian Hamiltonian H (6) is standard and appears in numerous 2D topological systems. However, the symmetry of H in momentum space is non-trivial and leads to the formation of bulk Fermi arcs.
Symmetry: Parity and Hermiticity. The hydrodynamic operator H hyd in (6) is Hermitian, a sum of products of Hermitian Pauli matrices and real numbers. Likewise, the elastic operator H elas is skew-Hermitian (i.e., a Hermitian operator multiplied by i), Note that the hydrodynamic forces do not conserve momentum, but the effective hydrodynamic operator is Hermitian. This is because, in the low-Reynolds regime, the sum of the forces is proportional to the velocity (the drag force in (4)). Mathematically speaking, the imaginary unit factors of the time derivative (iω) and the spatial derivative (ik) cancel each other. For the same reason, the effective elastic operator is skew-Hermitian, reflecting the overdamped nature of elastic modes in this regime.
As for parity symmetry in k-space, the hydrodynamic part is odd and the elastic part is even, This follows from the parity of the interactions (1,2): Ω x and Ω y are odd functions of k whereas ω 1 , ω x and ω y are even (see App. A 3). It is important to note that unlike "odd" systems [35,38], H does not introduce any nonreciprocal coefficients of elasticity or viscosity. Here, non-Hermiticity is simply the outcome of hydrodynamics.
The spectrum: Dirac cones and Fermi arcs. The interplay of odd, Hermitian hydrodynamics and even, skew-Hermitian elasticity brings about distinctive topological signatures (Fig. 2). The spectrum of equations (5,6) exhibits two eigenfrequency bands, where , and all frequencies are measured in units of τ −1 hyd . The corresponding polarization eigenstates are (App. A 4) Without elastic forces ( Fig. 2A), a purely hydrodynamic system ( = 0) exhibits real spectrum of propagating phonon-like waves [22,32], ω + = −ω − = (Ω 2 x + Ω 2 y ) 1/2 . On the edge of the Brillouin zone there are six Dirac points where the hydrodynamic interaction vanishes, Ω x = Ω y = 0. At a Dirac point, negative and positive bands kiss, ω + = ω − = 0, forming a graphenelike double cone [42] (App. A 6). The introduction of skew-Hermitian elasticity at non-zero breaks the symmetry. Fig. 2B shows the real part of the bands, Re(ω + ) = − Re(ω − ), for a hydroelastic number = π (the whole spectrum is shown in Fig. 7). Six bulk Fermi arcs-four S-shaped and two vertical lines-emerge from the Dirac points. Along the arcs, the real parts of the bands merge at Re(ω + ) = Re(ω − ) = 0. Each arc is therefore an openended isofrequency contour joining a pair of isolated exceptional points (ExPs). As shown below, these points are topological defects of opposite ± 1 2 charges (App. A 7).
Singularities and bifurcations. At the ExPs, the spectral bands (7) and their corresponding eigenstates (8) simultaneously coalesce (Fig. 2B): the bands are purely imaginary, ω + = ω − = −i ω 1 , and the coalescing eigenstates, e + k = e − k , are either right-or left-circular polarizations, reflecting the chirality of the topological charges (see App. A 7). ExPs occur where the determinant in (7) vanishes, (ν + = 0 or ν − = 0) at At these double branching singularities, of both Re(ω) and Im(ω), the spectrum becomes gapless. In contrast to the Dirac points whose eigenspaces are two-dimensional, the eigenstates at the ExPs are parallel, signifying a reduction of the eigenspace dimension to one. Tuning the hydroelastic number advances the ExPs along 1D trajectories, from the Dirac points at = 0 to the corners or the center of the Brillouin zone at = ∞ (App. A 7). Along these trajectories, the bands exhibit square-root singularities (Fig. 3), where s ≡ Ω x /ω y = −Ω x /ω y is a 1D coordinate defined by (9). The Fermi arc is the branch-cut of the square root (10) stretching between the bifurcation transitions at the ExPs, s = ± . The square-root singularity reflects strong level repulsion between the bands, compared to the linear opening of the gap at the Dirac cone [9].
Bulk Fermi Arcs are generic.
As mentioned, the bulk Fermi arcs are topological signature of non-Hermiticity [7] in the driven elastic lattice-unrelated to the more common surface Fermi arcs induced by Weyl points in 3D Hermitian systems [24][25][26][27]. The bulk arcs are a direct outcome of the broken symmetry, at = 0, when hydrodynamics and elasticity with opposite symmetries are mixed. Thus, the bulk arcs are a generic feature of the system, which are indeed observed in other lattice symmetries ( Fig. 9) and for all directions of the driving flow. This is formally shown by a small-expansion around the Dirac cone of a purely hydrodynamic system (App. A 8), demonstrating that solutions to the ExP double condition (9) always exist. Thus, the pair of ExPs and the bulk arc that stretches between them are generic.
Topological charges. The ExPs are topological defects: a closed eigenfrequency loop encircling an ExP cannot shrink to a point without passing through the ExP. The corresponding topological charge can be found from the vorticity of the band-gap [4], ∆ω = ω + − ω − (Fig 4), Since the band-gap is ∆ω = 2 √ ν + ν − (7), the vorticity with a sign corresponding to the left-or right-handed chiralities of the gradient field (App. A 10). The vorticity is determined by the branch-cuts, where the jumps of arg ν ± by ±2π yield the topological charges q ± = ± 1 2 . The charges and their opposite chiralities originate from the square-root singularity of the Riemann surface (10), and reflect the spinor-like nature of the polarization eigenstates, which accumulate a ±π phase when circling the ExP and passing through the branch-cut of the arcs.
Likewise, the charges can be computed from the integral of the Berry phase [43] (see App. A 9). Berry's connections are the vectors and Berry's phases γ ± are the loop integrals The corresponding charges, q ± = γ ± /(2π) = ± 1 2 , are determined by the jump of the phase (13) at the branch-cut, as in the case of the vorticity (11).
The density of states and its singularities. Projecting the Riemann surfaces of the spectral bands ω ± onto the complex frequency-plane reveals the density of states g(ω) (Fig. 5), The Dirac points (green) split between the two banks of the branch cut. Notable are sharp-edged ridges of the density merging at logarithmically diverging summits, akin to van Hove singularities in Hermitian systems [44]. Level repulsion shows in low density of states around the ExP singularities in Fig. 5. The analytic density of states g(ω) is similar to the one obtained from spectra of simulated lattices, with deviations owing to finite size and boundary conditions (Fig. 10).

DISCUSSION
Realizations and physical limitations. Four physical conditions determine the experimentally accessible regime: I. The topology is most notable when the elastic and hydrodynamic forces are comparable, around ∼ 1. Using (3), this condition can be written as κ ∼ ηu( /a) 3 (with Stokes law γ ∼ η ).
III. To avoid lattice "melting", the thermal fluctuations should be kept small. A Lindemann-like criterion implies that the fluctuations are smaller than the particle size, where k B T is the thermal energy scale and the fluctuations of the 2D elastic lattice increase logarithmically with size L [45]. This condition is recast as a lower bound, κ ≥ (k B T / 2 ) ln(L/a).
Examination of the bounds in the −κ plane (Fig. 6) shows that the regime they define is easily accessible in standard soft matter settings. The limits are drawn for the parameters of a prototypical microfluidic system [31].
We see that all conditions I-IV are met by the optimal spring constant κ ∼ 10 2 −10 3 k B T/µm 2 (from condition I), corresponding to interaction energies of order 10 2 −10 4 k B T between micron-sized particles. Physical settings of microscopic bodies whose interactions are of order ∼10 3 k B T are characteristic to various soft matter systems, in particular colloids [45]. For example, spring constants of about ∼10 2 −10 3 k B T/µm 2 were measured in colloidal gels [46], small lattices [47] and 2D colloidal crystals [48][49][50][51]. Note that elasticity can be replaced by any other conservative central force, such as magnetic [52] and electrostatic [53] interactionsall retain the Hermiticity and parity symmetries of elasticity and enter the equations of motion as an effective spring constant. Other relevant systems include lattices of DNA-coated colloids [54,55], colloidal lattices embedded in gels, and self-assembling microfluidic crystals [56,57].
To observe the exceptional topology, such soft matter crystals should be subjected to viscous flow, for example in a quasi-2D microfluidic channel. The lattices can be driven by the surrounding flow, slowed down by friction [31] or held in place by boundary conditions or external forces, for example, by dynamic holographic optical tweezers [58]. What matters is the velocity of the lattice relative to the embedding fluid, which induces the hydrodynamic interactions. The geometry of the lattice can be easily adapted. In our model system (Fig. 1A), we assumed for simplicity that the particles are joined by thin struts. However, the struts can be replaced by other particles, forming for example kagome, honeycomb or square-kagome lattices. This would augment the non-Hermitian spectrum with additional optical branches.
Moreover, long-range magnetic or electrostatic interactions obviate the need for direct mechanical connection among the particles. So one can leave out the struts and form, for instance, large 2D lattices of magnetic colloids [59], disks [60], hydrogel particles [61], or droplets [62], with effective springs induced by magnetic dipole-dipole interactions [52]. The steady state configuration of such a lattice is governed by the confining boundaries, which can be crafted to allow the passage of viscous flow.
Outlook and summary. The proposed non-Hermitian model system can also be realized also in the sub-micron regime. Shrinking the system's dimensions (Fig. 6) tightens the thermal stability bound on κ (condition III), which approaches the optimal value (condition I) at a particle size of about ∼ 100 nm. But the bounds may be relaxed by increasing ηu (shifting up the dashed line), for example using higher pressures in hardened microfluidic chambers [23,30]. Tentative candidate systems in this regime include lattices made of DNA [55,63], protein [64][65][66] and patchy nano-particles [49,67,68]. At the other extreme, the system can be expanded up to the sub-millimetric regime, as long as the effective springs remain soft. Minding the omnipresence of both periodic lattices [69,70] and low-Reynolds flows [71][72][73] in biology and biotechnology, one may speculate that non-Hermitian topological phenomena are also widespread in these realms.
To sum, the reported skewed topology with bulk Fermi arcs ending at isolated exceptional points was recently observed in a photonic crystal with radiation loss [7]. Here, we proposed a blueprint a controllable realization of this exotic topology in ordinary, passive soft matter, where the analogue of radiation loss is viscoelastic dissipation. The performed analyses showed the ubiquity of exceptional points and bulk Fermi arcs in generic elastic lattices driven by viscous flow. This suggests that non-Hermitian topology-which is robust against defects and noise-can be used to develop devices in the low-Reynolds regime. While these directions are intriguing, our main objective here was to demonstrate that simple viscoelastic systems can form an easily-accessible playground to investigate fundamental features of topological matter in the overdamped low-Reynolds regime, typical to soft and living matter. Future directions include examination of skin modes [6] and 3D geometries, where the long-range hydrodynamic interaction is much stronger than in 2D. In addition, one could explore composite architectures combining passive and active soft matter, also at the nanoscale [74,75].

ACKNOWLEDGMENTS
This work was supported by the Institute for Basic Science, Project Code IBS-R020. I thank Tsevi Beatus, Roy Bar-Ziv, Sam Safran, Hyuk Kyu Pak, Vincenzo Vitelli, William Irvine and their coworkers for essential comments.
Appendix A: Derivation of the model

Dynamics of hydroelastic lattices
The following is a brief description of the derivation. An elastic lattice is moving in the x-y plan of a thin 2D fluid layer at a velocity u relative to the fluid. In this quasi-2D geometry, the narrow dimension is z (perpendicular to the page in Fig. 1). The lattice is made of particles of size joined by thin elastic rods of average length a and spring constant κ. The viscous drag on each particle is γu, where γ is the friction coefficient (inverse of the mobility). Therefore, this non-equilibrium steadystate of uniformly moving lattice requires driving forces F (force per particle) that inject momentum and energy to compensate for the dissipative friction forces.
The particles' motion with respect to the fluid induces dipolar perturbations with a velocity field decaying as the inverse square of the distance, ∼u ( /r) 2 , and this dipolar flow fields give rise to collective hydrodynamic interactions [22,[30][31][32]. The hydrodynamic force f hyd ij exerted by the j th particle on the i th particle is (in x, y components), where the positions of the dipoles are R i = (X i , Y i ), and R ij = R i − R j = (X ij , Y ij ) are the distance vectors. In polar coordinates, R ij = (R ij , θ ij ), where R ij = |R ij | and θ ij the angle. Equation (1) is a compact form of (A1).
The coupling constant Λ scales as the strength of the dipoles, Λ ∼ u 2 , where is the size of the particle, and u its velocity relative to the fluid. This defines the typical timescale τ hyd of the hydrodynamic interaction, the time it takes a perturbation to traverse a distance a at a sound velocity c s ∼ Λ/a 2 (in the continuum long wavelength limit). The timescale τ hyd depends on the physical forces driving the hydrodynamic interaction. For example, in a lattice of particles sedimenting in a quasi-2D fluid, the relative velocity scales as u ∼ ∆ρ g 2 /η, where ∆ρ is the density difference and η the fluid's viscosity.
In the quasi-2D flow of squeezed droplets, the coupling is Λ = 2 Ku, where K = u particle /u fluid ≤ 1 is the ratio of the velocities of particle and the surrounding fluid that drags it in the channel [22,31]. The elastic Hookean forces are described by harmonic springs with constant k and equilibrium lengthR ij , the lattice constant a. The elastic forces are proportional to the change of the length, ∆R ij = R ij −R ij , where n ij is a unit vector in the direction of the distancē R ij . Equation (2) is (A3) expressed in polar coordinates.
In the overdamped regime, the elastic relaxation time is This is the typical decay time of vibrational modes in the absence of hydrodynamic driving force. The hydroelastic number measures the relative significance of elastic and hydrodynamic interactions (3), Unlike odd viscosity [38] or odd elasticity [35], both hydrodynamic and elastic forces are non-circular, ∇ Rij × f ij = 0, as they conserve angular momentum.

Linear Expansion
In the overdamped low-Reynolds regime, the equations of motion are (4). At steady-state, the lattice interactions vanish by symmetry, j =i (f hyd ij + f elas ij ) = 0, and the lattice moves uniformly at a velocity u = F/γ. Expansion of the equations of motion (4) in small deviations of the lattice positions around the steady-state positions, r j = R j −R j , yields a linear dynamic equation, where r is 2N -vector of the N particle deviations, r i . The tensor H = H hyd + H elas , which is analogous to a Hamiltonian, combines contributions from hydrodynamics and elasticity. H is a 2N × 2N -matrix, composed of 2 × 2 blocks H ij that account for interactions between the i th and j th particles, where the hydrodynamic term is and the elastic one which are conveniently expressed in terms of Pauli's matrices σ x , σ y , σ z (and 1 the unity matrix), H ij elas = 1 2 (1 + sin 2θ ij σ x + cos 2θ ij σ z ) .
The diagonal terms ensure zero sums, H ii = − j =i H ij . Since the angles obey θ ji = π + θ ij , it follows from (A9) that H ji hyd is odd with respect to particle exchange (i↔j), while H ij elas is even, The mutual forces between the i th and the j th particles are f ij = H ij (r j − r i ) and f ji = H ji (r i − r j ). Thus, we verify that the hydrodynamic forces violate Newton's third law of momentum conservation, whereas the lastic forces obey it (see Fig. 1B), While the microscopic molecular forces in the fluid obey Newton's law, the hydrodynamic interactions are effective macroscopic forces that do not conserve momentum. The momentum is leaking through the walls and is compensated by the driving force (e.g., gravitation or pressure gradient).

Momentum space
To exploit the crystal symmetry, one represents the dynamics in the momentum space of the wave-vectors k = (k x , k y ) by Fourier transform. The linearized dynamical equations (A5,A9) are expanded in plane waves, such that the deviation of each particle from its mechanical equilibrium positionR j is r j (t) = e k e ik·Rj = e k e i(k·Rj −ωt) , where e k (t) = e k exp(−iωt) is a 2D polarization vector in k-space.

The equation of motion of e k (t) is a non-unitary
Schrödinger-like equation, whose eigenvector is e k with eigenfrequency ω (5), The "Hamiltonian" H is a 2×2-matrix expressed in terms of Pauli's matrices as The contributions of the long-range hydrodynamic to H (A12) are Fourier sums, whereR j =R j (cos θ j , sin θ j ) are the distances of the steady-state lattice positions from an arbitrary origin particle 0. Owing to the crystal's parity symmetry, we can rearrange the summation to be over pairs at inverse positions, ±R j , demonstrating that Ω x and Ω y are always real in a crystal, The elastic contributions in (A12) are sums over the neighbors connected by springs to particle 0, Again, owing to the crystal parity symmetry, the sums can be rearranged, demonstrating that ω x , ω y and ω 1 are all real, The hydrodynamic and elastic interaction in a triangular lattice are shown in Fig. 8.

Spectra
The eigenfrequencies ω are found by solving the secular equation corresponding to (A11,A12). There are two eigenfrequency bands, where ν ± are defined as In the hydrodynamics-dominated regime, 1, the spectrum (A15) is purely real, while in elasticity-dominated regime 1 the spectrum is purely imaginary, The polarization eignevectors, e + k and e − k , are where the eigenvalues ω ± are given in (A15). When normalized, the eigenvectors become

Circular polarization basis
One can represent the Hamiltonian in the basis of left and right circularly polarized unit vectors (which are the coalescing eigenstates at the ExPs), In this basis, the Hamiltonian (A12) becomes equation (6), The eigenvectors in this representation take the simple form
The expansion of (A13) around the Dirac point is linear where the gradients are The α j 1 and α j 2 in (A21) are the indices of the lattice points,R j = α j 1 a 1 + α j 2 a 2 , with the basis vectors, a 1 and a 2 . It is straightforward to verify that the gradients at the Dirac point are orthogonal, resulting in an elliptic cone ( Fig. 2A). It is also useful to note that, at the Dirac points, the elastic interactions, ω x , ω y and ω 1 , have extrema (minima, maxima, and saddles), and are therefore constant to first order in dk (Fig. 8). 9. The spectrum of a square lattice. The real (top) and imaginary (bottom) components of the eigenfrequencies, ω+ (left) and ω− (right), of a square lattice at = π. The black square is the boundary of the Brillouin zone. Shown are the Dirac points (green) and the ExPs with their + 1 2 (red) and − 1 2 (blue) charges. Bulk Fermi arcs are grey curves.

Exceptional points and Fermi arcs
Exceptional points (ExPs) occur when both eigenfrequencies (A15) and eigenvectors (A19) coincide, where the ± signs correspond to the paired ExPs. In the circular polarization basis, the eignevectors at the ExPs are the unit vectors, ExPs are positioned exactly where the determinant, ν + ν − , in (A15) vanishes. In other words, when ν + = 0 or ν − = 0, which amounts to the double condition (9), As one varies , the ExPs move along the 1D curve from the Dirac cones at = 0 to the Brillouin zone corners (or its center) at = ∞. The stretch of the trajectory between the two ExPs is the Fermi arc. Along these trajectories, the eigenfrequencies are where s ≡ Ω x /ω y = −Ω x /ω y is the coordinate along the trajectory, from s = 0, the Dirac point, through s = , the ExP, to the corner of the Brillouin zone s = ∞ (the four S-shaped arcs), or its center (the two straight arcs). The exceptional point is the position of a bifurcation transition where the branch-cut of the square root stretches along the Fermi arc.

Bulk Fermi arcs are generic
The bulk Fermi arcs observed in the triangular (Figs. 2,7) or the square (Fig. 9) lattices are a generic phenomenon of the hydroelastic lattice. They emerge as soon as symmetry is broken, at = 0, when the skew-Hermitian elasticity is introduced. To see this formally, one may examine the emergence of arcs for small . Then, the exceptional point condition (9) can be linearly expended in small deviations from the Dirac point, dk = k − k D (A21). The double condition (9) becomes a linear equation in dk, where the ± signs correspond to the two ExPs. Since the gradients are linearly-independent (A22), solutions always exist, and the positions of the ExPs are Hence, the arcs and the ExPs are generic. The bulk Fermi arc stretches between the two solutions, and its size grows linearly with .

Berry's phase and the topological charges
In the momentum representation, Berry's connections are the vectors [43], corresponding to the eigenvectors e ± k . For convenience, we write the normalized eigenvectors as Jones polarization vectors, In this notation, we find that the Berry connections are simply For the eignevectors (A18), the Jones parameters are c = √ ν + + √ ν − 2 2 (|ν + | + |ν − |) , The Berry phases γ ± are the path integrals To find the Berry phase of the ExPs, we integrate over a small circular counterclockwise path around each point. For example, consider the point where ν + = 0. In the vicinity of an ExP, c 1 2 , and the connections are therefore, Along the circular integration path, the phase arg ν + passes all quadrants. In particular, it passes through the branch cut, arg ν + = ±π, where it jumps by ±2π, depending on the gradients of ν + . Thus, we find that Berry's phase is γ ± = ±π. Likewise, at the other ExP, Berry's phase takes the opposite sign, γ ± = ∓π. Altogether, we find that the corresponding topological charges are q ± = γ ± /(2π) = ± 1 2 .

Vorticity and Topological Charges
Another way to find the topological charges is through the vorticity of the eigenfrequency band-gap, ∆ω = ω + − ω − , defined as [4] where ω ± are the given in (A13). Since, ∆ω = ω + −ω − = 2 √ ν + ν − , we find that the vorticity is As in the case of the Berry phase, the vorticity around the ExPs is determined by their branch-cuts. Namely, the jumps of arg ν ± by ±2π when passing through the branch cuts yield the topological charges q ± = ± 1 2 . and simulation (B), for = π. The simulated system is a 121 × 121 triangular lattice with periodic boundary conditions (a torus). The spectrum is obtained by Fourier analysis of the motion computed from the dynamical equation (4). The simulated g(ω) exhibits general similarity to the analytic solution. The incommensurately of the four-fold symmetry of the torus and the six-fold symmetry of the lattice results (C-D): g(ω) for square lattice, analytic (A) and simulation of a 121 × 121 lattice (B), for = π. The correspondence between the simulation and the analytic solution is excellent, owing to the commensurable four-fold symmetry of the lattice and the periodic boundary conditions, with some coarseness due to the finite size of the simulated system.