Parity-time phase transition in photonic crystals with \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$C_{6v}$$\end{document}C6v symmetry

We investigate the parity-time (PT) phase transition in photonic crystals with \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$C_{6v}$$\end{document}C6v symmetry, with balanced gain and loss on dielectric rods in the triangular lattice. A two-level non-Hermitian model that incorporates the gain and loss in the tight-binding approximation was employed to describe the dispersion of the PT symmetric system. In the unbroken PT phase, the double Dirac cone feature associated with the \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$C_{6v}$$\end{document}C6v symmetry is preserved, with a frequency shift of second order due to the presence of gain and loss. The helical edge states with real eigenfrequencies can exist in the common band gap for two topologically distinct lattices. In the broken PT phase, the non-Hermitian perturbation deforms the dispersion by merging the frequency bands into complex conjugate pairs and forming the exceptional contours that feature the PT phase transition. In this situation, the band gap closes and the edge states are mixed with the bulk states.


Results
two-level non-Hermitian model. Consider a two-dimensional photonic crystal consisting of equilateral triangular rods arranged as hexagonal clusters in the triangular lattice, as schematically shown in Fig. 1a. The wave equation in terms of the time-harmonic magnetic field H(r) (with the time convention e −iωt ) is given by where ε(r) is the complex dielectric function. Assume that the lattice possesses the parity-time (PT) symmetry, that is, the Hamiltonian commutes with the combined PT operator: (P T )L = L (P T ) , where P and T are the parity-inversion and TR operators, respectively. A necessary condition for the PT symmetry to hold is www.nature.com/scientificreports/ that the dielectric function satisfies ε(r) = ε * (−r) , where ε(r) ≡ ε r (r) + iε i (r) with both ε r and ε i being real. For a triangular lattice composed of dielectric rods, the PT symmetry is established by an antisymmetric arrangement of gain and loss on the rods in the unit cell (cf. Fig 1b). Denote γ and g(r) the strength and distribution of gain and loss, respectively, that is, ε i (r) ≡ γ g(r) , where γ ≡ max r |ε i (r)| and −1 ≤ g(r) ≤ 1 . The gain and loss are applied such that g(r) = −g(−r) and the convention γ ≥ 0 is adopted. In the present problem, we propose that ε i (r) changes sign by a rotation of 60 • about the unit cell center. If γ is zero, L is a Hermitian operator and ω is always real 47 .
Suppose that γ is small compared to the maximum of |ε r (r)| , and iγ g(r) is treated as a small perturbation to the dielectric function 38 . Using the binomial expansion for 1 ε(r) , the wave equation is approximated as In a two-level non-Hermitian model 46 , the operator L is separated into a Hermitian part L h ≡ ∇ × 1 ε r (r) ∇× and a skew-Hermitian part L s ≡ −∇ × iγ g(r) ε 2 r (r) ∇× . For a structure with C 6v point group symmetry, there exist doubly degenerate E 1 and E 2 states, with the polynomial representations {x, y} and {2xy, x 2 − y 2 } , respectively 48 . In a triangular lattice with C 6v symmetry, the eigenmodes of the E 1 and E 2 symmetries are also identified at the Brilloin zone center (the Ŵ point) 21 , which are referred to as the p and d orbitals, respectively 13 .
Effective Hamiltonian. Based on the tight-binding approximation for the triangular lattice, the eigenfield around the p and d bands at the wave vector k near the Ŵ point is expressed as 21 where H (1) and H (2) ( H (3) and H (4) ) are the normalized magnetic fields of the E 1 ( E 2 ) state for a single unit structure with C 6v symmetry and α j ( j = 1, 2, 3, 4 ) is the weighting coefficient. Here, m = 0 denotes the center cell, m = 1, 2, . . . , 6 its nearest neighboring cells with the centers at r m , and e ik·r m is included to satisfy the Bloch theorem: H k (r + R) = e ik·R H k (r) on the lattice. The states H (1) , H (2) , H (3) , and H (4) form a basis for the four-band subsystem, which is denoted as p x , p y , d 2xy , d x 2 −y 2 . Using Eq. where H is the Hamiltonian of the four-band system and I is the 4 × 4 identity matrix.
In the present two-level model, the effective Hamiltonian for the Hermitian operator L h on the basis p x , p y , d 2xy , d x 2 −y 2 in the linear order of k is given by (see Methods C.1) where A is a real quantity 13 . The Hermitian Hamiltonian H h has the same form as in the photonic structure with C 6v symmetry 26,49 . The effective Hamiltonian for the skew-Hermitian operator L s on the same basis is given by (see Methods C.2) where N 1 and N 2 are real quantities. The effective Hamiltonian for the non-Hermitian operator L , given by H = H h + H s , has a block form of the PT symmetric matrix with a generic (2, 2) parity operator 50 . A similar form of the non-Hermitian Hamiltonian was also present for the photonic structures with C 4v symmetry 38 . pt phase transition. Based on Eqs. (5), (6), and (7), the eigenfrequency of the four-band non-Hermitian system can be solved to give the following formula:  21 . In general, ω p = ω d and a gap between the p and d bands is opened. For γ > 0 , the eigenfrequencies of the four bands may not be real. Denote the discriminant of the square root in Eq. (8) as If γ is small enough such that � > 0 , the eigenfrequencies are real and the PT symmetry is unbroken. On the other hand, if γ exceeds a threshold value such that � < 0 , the p and d bands are merged into complex conjugate pairs and the PT symmetry is broken. In this situation, the real parts of the p and d bands overlap at the their quadratic mean frequency: ω 0 = ω 2 p + ω 2 d /2 and the gap in between is closed.
The exceptional point that features the transition between an unbroken and a broken PT phase is determined by the quadratic equation: = 0 . A collection of exceptional points that separate the regions of real eigenvalues and complex conjugate pairs form a contour in the wave vector space, the so-called exceptional contour 38 or exceptional ring 35 . In the present problem, the PT phase transition is considered a balance between the strength of gain and loss ( γ ) and the deviation from the four-fold degeneracy either in frequency ( ω 2 p − ω 2 d ) or in wave number ( k 2 x + k 2 y ). In case the accidental degeneracy at the Ŵ point is attained, where ω p = ω d and k x = k y = 0 , we have � = −4γ 2 N 2 1 + N 2 2 < 0 for any nonzero γ and the eigenfrequencies of the p and d bands are always complex conjugate pairs. In this situation, the PT symmetry is spontaneously broken in the presence of an infinitesimal amount of gain and loss 34 . The thresholdless PT symmetry breaking also occurs at the degenerate point with two-fold [36][37][38] or three-fold 35 degeneracy. It can be shown that the thresholdless PT symmetry breaking at the degenerate point is valid for a non-Hermitian system with finite decoupled bands (see Methods D).
If the non-Hermitian system is in the unbroken PT phase, that is, the eigenvalues are real, the perturbed frequency due to the presence of balanced gain and loss can be estimated by the perturbation method. For this purpose, a small number ǫ is defined as ǫ ≡ γ /ε m , where ε m ≡ max r |ε r (r)| . The frequency shift of the perturbed system is given by (see Methods E) where H 0 is the magnetic field at the eigenfrequency ω 0 for the unperturbed system (without gain and loss). Note that the first term on the right-hand side of Eq. (9) is purely imaginary provided that ω 0 is real (which is true when the unperturbed system is Hermitian). In the unbroken PT phase, the eigenvalues are real, which dictates that the first-order correction to the eigenfrequency is zero and the frequency shift due to the presence of balanced gain and loss is of second order: �ω = O ǫ 2 . A similar form of the frequency shift can be obtained when the perturbation method is applied on the electric field: �ω = − ω 0 where E 0 is the electric field with the eigenfrequency ω 0 51 .

Discussion
Broken pt phase. Figure 2 shows the band structure of TM modes along the high symmetry points of the first Brillouin zone for the PT symmetric photonic crystal with a 0 /R = 3 . Here, the TM modes are characterized by the in-plane magnetic field components ( H x and H y ) and the out-of-plane electric field component ( E z ). If the gain and loss is not present ( γ = 0 ), the photonic crystal possesses C 6v symmetry 21 . In this situation, the photonic crystal has two equivalent lattice descriptions (as a triangular lattice or a honeycomb lattice), which corresponds to the transition point between a trivial and a nontrivial topological phase 13 . For a small strength of gain and loss ( ε d = 12 + 0.1i ), the eigenfrequencies deviate from the lossless case up to O ǫ 2 [cf. Eq. (9)]. The real parts of the second to fifth bands (p and d bands) are nearly degenerate at the Ŵ point (Fig. 2a), where the band dispersion is mainly described by a double Dirac cone [cf. Eq. (8)]. In the presence of balanced gain and loss, however small, the non-Hermitian system is spontaneously broken at the degenerate point 34 . The thresholdless PT symmetry breaking is manifest on the emergence of tiny imaginary parts of the p and d bands around the Ŵ point (Fig. 2b).
The band structure of the PT symmetric photonic crystal at a 0 /R = 2.75 for a larger strength of gain and loss ( ε d = 12 + 1.8i ) is shown in Fig. 3. The p and d bands merge into complex conjugate pairs in the region near the Ŵ point, where the real parts of the eigenfrequencies overlap at their quadratic mean frequency (Fig. 3a). Meanwhile, the imaginary parts of the degenerate bands tend to split from the Ŵ point to give two loops in the overlapping region (Fig. 3b). In this situation, the Dirac dispersion of the non-Hermitian system is deformed into flat bands near the center (Fig. 3c). The exceptional points that separate the regions of real eigenvalues (inside) and complex conjugate pairs (outside) form two contours in the wave vector space (Fig. 3d), which are characterized by the quadratic equation: = 0 [cf. Eq. (8)]. Near the Ŵ point, the dispersion surface is approximately isotropic with respect to k x and k y , and the exceptional contours are roughly circles.  Figure 4 shows the band structures and eigenfields for the PT symmetric photonic crystal with ε d = 12 + 0.6i , where the bulk states are in the unbroken PT phase with real eigenvalues. For a 0 /R > 3 , the p bands lie below the d bands ( ω p < ω d ) (Fig. 4a), that is, the dipole modes with odd parity (Fig. 4d,e) have the eigenfrequency lower than the quadrupole modes with even parity (Fig. 4b,c), and the photonic crystal has a trivial topological phase as in the Hermitian case. For a 0 /R < 3 , on the other hand, the p bands lie above the d bands ( ω p > ω d ) (Fig. 4f), and the dipole modes (Fig. 4g,h) have the eigenfrequency  www.nature.com/scientificreports/ higher than the quadrupole modes ( Fig. 4i,j). In this situation, the band inversion (or parity inversion) between the p and d bands occurs, and the photonic crystal has a nontrivial topological phase as in the Hermitian case.
edge states. Figure 5 shows the dispersions and wave propagations of the edge states at the interface between two PT symmetric photonic crystals, with a trivial topological phase ( a 0 /R = 3.35 ) above the interface and a nontrivial phase ( a 0 /R = 2.75 ) below. In Fig. 5a, both crystals are in the unbroken PT phase ( ε d = 12 + 0.6i ). For a small strength of gain and loss, the edge states, which lie inside the common band gap for the two crystals, can exhibit an entirely real spectrum along the k x direction. Notice that the edge states do not cross each other at the Ŵ point and a small gap is opened in between, which means that the edge states may not be gapless 26 . This feature is attributed to the breaking of C 6v symmetry at the boundary between two photonic crystals with different a 0 /R , which is analogous to the effect of magnetic impurity that breaks the TR symmetry in a QSH system and opens a gap in the edge state 52 . As there is no point degeneracy in the edge states, the thresholdless PT symmetry breaking as in the bulk states does not occur. In Fig. 5b, a circularly polarized magnetic dipole composed of two perpendicular in-plane dipoles with 90 • phase difference is placed at the interface between a trivial phase (above the interface) and a nontrivial phase (below the interface) to excite the surface wave at a frequency inside the common gap (cf. blue dashed line in Fig. 5a), where the field is evanescent on either side of the interface. Excited by the dipole with H 0 e −iωt x + iŷ , the right-handed wave propagates unidirectionally toward the right, which is consistent with the surface band dispersion with a positive k x (cf. blue dot in Fig. 5a). In particular, the surface wave is able to bend around sharp corners without backscattering. As the strength of balanced gain and loss increases, the gap in the edge states is reduced. At a certain value of γ , the edge states can be gapless in the real part, with the emergence of imaginary part at the crossing point. For a larger value of γ in the nontrivial phase ( ε d = 12 + 1.2i ), the edge states merge into a complex conjugate pair, while the bulk band gap is still open, as shown in Fig. 5c. The edge states exhibit flat bands in the region near the Ŵ point in a similar manner as the broken PT phase of the bulk states (cf. Fig. 3a). In Fig. 5(d), the left-handed wave is excited by H 0 e −iωt x − iŷ at the same frequency (cf. blue dashed line in Fig. 5c) where the edge state has a real eigenfrequency (although the spectrum is not entirely real). The edge state propagates unidirectionally www.nature.com/scientificreports/ toward the left, which is consistent with the surface band dispersion with a negative k x (cf. blue dot in Fig. 5c). As in Fig. 5b, the surface wave is able to bend around sharp corners without backscattering. This feature holds if the larger value of γ is given instead in the trivial phase. In Fig. 5a-d, the bulk states of both crystals are in the unbroken PT phase. The spin-polarized edge states retain the helical nature as in the Hermitian case, where the two states with opposite helicity (handedness) counterpropagate at a given edge and are robust against disorder.
For an even larger value of γ in the nontrivial phase ( ε d = 12 + 2i ), where the PT symmetry is broken, the merged bulk bands occupy the original band gap region and the edge states are no longer separated from the bulk states, as shown in Fig. 5e. In this situation, the dipole source at the same frequency (cf. blue dashed line in Fig. 5e) excites not only the surface state (still having a real eigenvalue at the exciting frequency) but also the bulk states. The surface waves are mixed with the bulk waves that are scattered over the whole region on the lower side, as shown in Fig. 5f. If, on the other hand, the PT symmetry is broken in the trivial phase, the surface waves are to be mixed with the bulk waves on the upper side.
In conclusion, we have investigated the PT phase transition in photonic crystals with C 6v symmetry, with balanced gain and loss on dielectric rods in the triangular lattice. The dispersion of the PT symmetric system was    where H is a 4 × 4 matrix with the entries H ij ( i, j = 1, 2, 3, 4 ) and I is the identity matrix.
B. Matrix representations of the E 1 and E 2 modes in C 6v symmetry group. For the structure with C 6v point group symmetry, there exist doubly degenerate E 1 and E 2 states with the polynomial representations {x, y} and {2xy, x 2 − y 2 }, respectively 48 . For the E 1 state, the matrix representations of the elements (used in the present study) in C 6v symmetry group (cf. Fig. 6b) are given by 53 For the E 2 state, the corresponding matrix representations are given by  www.nature.com/scientificreports/ As the magnetic fields of the basis states H (i) ( i = 1, 2, 3, 4 ) are highly localized in the single unit structure 24,54 , the field couplings between the neighboring cells will be much weaker than those in the same cell. where R H (i) (r) ≡ R H (i) R −1 r . In the present PT symmetric configuration, the real part of the dielectric function satisfies R ε r (r) = ε r (r) , while the imaginary part changes sign under the C 6 operation (rotation by 60 • ): Using Let P 0 be the matrix representation of the parity operator P on the basis p x , p y , d 2xy , d x 2 −y 2 . Based on the C 2 operation for the E 1 and E 2 states [cf. Eqs. (18) and (19)], we have which is a generic (2, 2) parity operator 50 . Since H s commutes with P 0 T (a consequence of the PT symmetry for both H and H h ), it can be shown that 50 and H s has the block form where A 2×2 , B 2×2 , C 2×2 , and D 2×2 are real matrices. This means that G (ij) 0 and G (ji) 0 are purely imaginary for i = 1, 2 and j = 3, 4 . Using the C 6 operation in Eqs. (18), (19), and (28)     www.nature.com/scientificreports/ D. thresholdless pt symmetry breaking. Suppose that n frequency bands for a PT symmetric non-Hermitian system can be decoupled from the rest of the system 38 and their eigenstates form a basis in the subsystem.
Proposition If the n bands in the decoupled subsystem are degenerate at a specific point, then the subsystem is Hermitian at that point.
Proof Let φ i ( i = 1, . . . , n ) be the independent eigenstates of the Hamiltonian H for the subsystem with the same eigenvalue 0 , that is, Hφ i = 0 φ i . Since any state in the subsystem can be written as a linear combination of the basis: i c i φ i , we have H = 0 I , where I is the n × n identity matrix. From H = P T P T H = P T HP T = P T HT P = P H * P = P * 0 I P = * 0 I , we know that 0 = * 0 and thus 0 is real and H is Hermitian.
By contraposition, if H is non-Hermitian then the eigenvalues of the PT symmetric subsystem cannot be all equal. At the n-fold degenerate point where all eigenvalues in the subsystem have the same real part, some eigenvalues should have nonzero imaginary parts (otherwise, they will be all equal). This means that the PT symmetry is broken at the degenerate point (without threshold) for the non-Hermitian system stated above. e. perturbation method. Let the presence of gain and loss in the non-Hermitian system be regarded as a perturbation to the Hermitian system 40,41 . By defining a small number ǫ ≡ γ /ε m , where ε m ≡ max r |ε r (r)| , the eigensystem is written as where L 0 = L h = ∇ × 1 ε r (r) ∇× and L 1 = 1 ǫ L s = −∇ × iε m g(r) ε 2 r (r) ∇× . Let H 0 be the eigenfield at the frequency ω 0 for the unperturbed system, that is, L 0 H 0 = ω 2 0 c 2 H 0 . The eigenfield for the perturbed system is expanded in a series as H = H 0 + ǫH 1 + ǫ 2 H 2 + · · · . Likewise, the eigenfrequency is expanded as ω = ω 0 + ǫω 1 + ǫ 2 ω 2 + · · · . Using these expansions in Eq. (39), we have The perturbed eigenfield H is normalized according to �H 0 |H� = �H 0 |H 0 � , so that the orthogonality condition �H 0 |H n � = 0 holds for n > 0 . The first-order term of the eigenfrequency is given by The first-order correction to the eigenfrequency is then obtained as