Broadband circular dichroism in chiral plasmonic woodpiles

The circular dichroism (CD) of a material is the difference in optical absorption under left- and right-circularly polarized illumination. It is crucial for a number of applications, from molecular sensing to the design of circularly polarized thermal light sources. The CD in natural materials is typically weak, leading to the exploitation of artificial chiral materials. Layered chiral woodpile structures are well known to boost chiro-optical effects when realized as a photonic crystal or an optical metamaterial. We here demonstrate that light scattering at a chiral plasmonic woodpile, which is structured on the order of the wavelength of the light, can be well understood by considering the fundamental evanescent Floquet states within the structure. In particular, we report a broadband circular polarization bandgap in the complex band structure of various plasmonic woodpiles that spans the optical transparency window of the atmosphere between 3 and 4 \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$\upmu$$\end{document}μm and leads to an average CD of up to 90% within this spectral range. Our findings could pave the way for an ultra-broadband circularly polarized thermal source.


Introduction
The manipulation of light through optical elements such as lenses, color filters, and polarizers goes back to ancient times and can ubiquitously be found in the living world [1,2]. The limits that naturally available materials impose on light manipulation can be overcome by nano-structuring matter in form of for example photonic crystals, dielectric geometries that are structured on the order of the wavelength of the light and lead to interference effects [3,4]. A different strategy to obtain new optical materials employs plasmonic metals structured on a deeply sub-wavelength scale, so-called metamaterials [5]. A vast number of chiral photonic crystal and metamaterial designs [6] have been suggested to yield custom-designed chiro-optical effects, such as strong circular dichroism [7][8][9][10][11][12], optical activity [13][14][15][16][17][18], and orbital angular momentum generation [19][20][21]. Next to a number of applications in sensing, catalysis, and chiral light generation, when realized on a micrometer length scale, these chiral geometries can be engineered to yield circularly polarized thermal emission by application of Kirchhoff's law [22][23][24][25][26][27][28][29][30][31][32].
A particularly promising and relatively easy to manufacture (down to sub-micrometer length scales) geometry is the chiral woodpile [34][35][36][37][38] as illustrated in Fig. 1. When realized as a classical deeply sub-wavelength metamaterial, the chiral woodpile, which is evidently not based on local meta-atoms, does not show a strong chiro-optical response due to the mismatch between the screw-axis pitch and the vacuum wavelength. A hexagonal chiral woodpile realized as a lossy photonic crystal with big index contrast (using a semiconductor with = 8.9 ), on the other hand, is predicted to produce circularly polarized thermal emission within a broad band (gap-to-midgap ratio Δ ∕⟨ ⟩ ≈ 1∕6 ) [22]. While a more sophisticated fully three-dimensional geometry with = 12 is predicted to improve the bandwidth to a gap-to-midgap ratio of 1/3 [12], we here instead consider a plasmonic hexagonal woodpile structured on the order of the wavelength that leads to a strongly circularly polarized thermal emission within a band with a gap-to-midgap ratio ≳ 2∕7 , focusing on the optical transparency window of the atmosphere between 3 and 4 μ m. A similar structure has been previously considered to engineer linearly polarized thermal emission [36]. Since the spatial periodicity of such a plasmonic crystal (PC) is on the order of the wavelength, standard metamaterial homogenization techniques cannot be applied. Looking at the field patterns can shed some light on the chiral plasmonic excitations, as for example done for a helix PC [39].
The simulated scattering of an outside source, however, always corresponds to an a-priory unknown superposition of self-consistent modal solutions within the PC. A physically inspired set of such modal solutions is provided through evanescent Floquet states [40][41][42]. These states can be thought of as a natural complex generalization of the Bloch modes discussed in the canonical band structure picture employed in the photonic crystal theory [3]. They thus allow to account for base materials with general optical permittivity and permeability and span the solution space of the monochromatic Maxwell equations within finite slab-like geometries with periodicity. While the calculation of evanescent Floquet modes is not implemented in standard simulation tools based on the finite element or finite difference methods, they can be naturally obtained through a plane-wave approach [40] for photonic crystals. For plasmonic materials, such an approach is, however, very inefficient as the generally poor convergence behavior of the plane-wave basis becomes a particularly big problem at metal-dielectric interfaces [41], where the permittivity changes its sign. In simple geometries, which are homogeneous in the propagation direction, such as the lamellar grating [43], the fishnet structure [42], or hyperbolic aligned wire media [44,45], the computation simplifies to an essentially one-or two-dimensional problem. More generally, in a geometry that is made by a number of slices that are individually homogeneous in the propagation direction, the evanescent Floquet modes are efficiently calculated as the eigen-solutions of the corresponding transfer matrix through the unit cell [46][47][48][49][50].
We here show that the fundamental Floquet modes responsible for propagating energy through (and absorbing it within) a single lamellar grating layer are sufficient to calculate the two fundamental Floquet modes of a hexagonal chiral woodpile made of a sequence of equal lamellar grating layers, copied through a C 6 2 screw rotation, as illustrated in Fig. 1d. These fundamental Floquet states of the chiral woodpile in turn describe the scattering physics well and reveal the origin of a broad circularly polarized stop band in the absorption spectrum, which is induced through a circularly polarized pseudo-bandgap. We use these two Floquet modes as a scattering basis within the fundamental Bragg order to efficiently optimize the broadband circular dichroism and derive general design principles based on the thickness of the dielectric region d 1 and the chiral pitch c only (Fig. 1). We show that the obtained absorption spectra are surprisingly accurate considering the strong approximations made and agree well with full-wave simulations. Our findings suggest experimental investigation of the designed  [33] is generated through an active C 3 (120 • ) rotation of the coordinate frame R after each lamellar layer. The unit cell transfer matrix is thus given by T UC = (RT) 3 , with the transfer matrix through one lamellar grating layer T woodpiles for chiral thermal emission. More fundamentally, they demonstrate the power of the concept of evanescent Floquet states to understand light propagation in finite slabs of artificially structured materials.

Methods
Evanescent Floquet states span the vector space of solutions of Maxwell's equations within a finite slab of a periodic material [40]. While in theory, a countably infinite number of such modes exists at each frequency, a small number usually suffices in practice, making them an invaluable tool to not only predict but most importantly understand the physical origin of a scattering experiment. If a periodic optical material contains base materials with loss, such as plasmonic metals, even the fundamental Floquet modes do not fit into the standard band structure picture with real-valued wave vectors [3]. Instead, a complex band structure picture has to be adopted [40,41]. Calculating these modes for general geometries remains challenging and is currently not possible with established software packages (based on for example finite differences or finite elements). We here present a semianalytical method to efficiently compute the complex band structure for the fundamental evanescent Floquet states of a hexagonal chiral woodpile. While our method cannot accurately predict scattering observables, it instead 1. reveals the main physical mechanism of scattering even for frequencies close to the first Wood anomaly substantially above the homogenization regime, and 2. allows for easily optimizing a desired behavior (broadband circular dichroism in this manuscript) within a large parameter space.

Single layer scattering
We start by calculating the Floquet states within each layer, which is homogeneous in z direction and whose local coordinate frame is defined in Fig. 1c. The procedure is based on a symmetry simplification of the well-known lamellar grating equation [43]. At normal incidence, the electric (TE) or magnetic (TM) field points in y-direction, and is thus antisymmetric with respect to the mirror operation y that maps y ↦ −y . The relevant scalar (monochromatic) field F(r) is the y-component of the electric (TE) or magnetic (TM) field. For the whole chiral slab, since we are interested in plane-wave excitation at normal incidence, we can characterize the field as odd with respect to the two-fold rotation symmetry around the z-axis at the center of the unit cell C 2z , that coincides with and results from multiple application of the C 6 2 -axis shown in cf. Fig. 1d. The scalar field F(r) is, therefore, symmetric under x = C 2z y at the center of both the dielectric and the metal domain, shown as dashed lines in Fig. 1c.
Due to the homogeneity of the single layer in z-direction, the Floquet modes are plane-wave-like. The general monochromatic solution of Maxwell's equations that satisfies the symmetry requirements in the domain = 1, 2 , where Here, c ∈ ℂ is a complex coefficient, q ∈ ℂ the wave number in the propagation direction, x the center of the domain, and k = ± √ k 2 0 − q 2 the lateral wave number given by the material dispersion relation, with k 0 ∶= ∕c 0 the vacuum wave number (with the angular frequency and c 0 the speed of light). The solution at the interface between the two domains additionally requires the tangential components of the electric and the magnetic field to be continuous, which leads to with the wave impedance Z = k ∕k 0 (TE) and Z = k ∕( k 0 ) (TM) and the (generally complex-valued) optical phase = k d ∕2 . A countably infinite number of Floquet solutions is thus obtained from solving the root equation This is evidently a transcendental problem as all and Z implicitly depend on q through the respective material dispersion relation. For the metal-dielectric structure and the wavelength range under consideration, however, a good guess can be obtained analytically. For this, let us first discuss the function (q) in the complex plane. It inherits two branch points from the two root functions of k at q = √ k 0 .
We now use the fact that these two roots are well separated (compared to k 0 ) in the mid-infrared, where ‖q 2 ‖ ≫ ‖k 0 ‖ for all metals. Due to the impedance mismatch between the two regions, the intensity of the low-order Floquet modes is either concentrated in the dielectric or in the metal region. In other words, solutions will be either found relatively close to q 1 or close to q 2 (and far away from the other branch point). Let us start with the first case, where the intensity is concentrated in the dielectric. Since the solution is far away from q 2 , the phase 2 has a large, positive imaginary part ℑ[ 2 ] ≫ 1 . The root function thus simplifies to ∶= Z 1 sin 1 cos 2 + Z 2 sin 2 cos 1 .
If the fields laterally concentrate in the metal domain, the phase 1 has a large, negative imaginary part ℑ[ 1 ] ≪ −1 , leading to Since ‖k 1 ‖ ≫ ‖k 2 ‖ , and ‖ 1 ‖ ≪ ‖ 2 ‖ , we obtain ‖Z 1 ‖ ≫ ‖Z 2 ‖ irrespective of polarization. The approximate sequence of roots for both polarizations is thus In the remainder of this paper, we only consider the fundamental TE/TM air mode with = = 0 , which is sufficient to qualitatively explain the observed scattering physics for small metal fill fractions 1 and a lateral lattice constant below the first Wood anomaly ( a < ). The successful application of this crude approximation can be understood by considering the following two facts: First, higher-order dielectric modes have little intensity in the fundamental Bragg scattering order (the DC Fourier component in x). On the other hand, the metal modes exhibit a weak coupling to a vacuum plane wave due to the strong impedance mismatch. The exact roots of the fundamental mode can be found using a standard Newton procedure with (4) as an initial estimate.
The fundamental and first-order modes of a platinum-air grating are shown in Fig. 2. Higher-order TE modes, where the Newton formalism does not converge with the analytical guess values, can be obtained using a global contourintegral method on a disk in the complex plane around the approximate higher-order root positions [51]. We illustrate the behavior of (q) , including the branch cuts, the position of the roots, and the contour integral method in Fig. 3.

Chiral woodpile unit cell transfer matrix and complex band structure
We conveniently calculate the transfer matrix within one unit cell using the single layer transfer matrix T and the rotation matrix R between neighboring layers, as illustrated in Fig. 1d. The layer transfer matrix connects the parallel components of the fields (E x , E y , H x , H y ) at the top of the layer to those at its bottom in its native coordinate frame, shown in Fig. 1c, d. We first introduce the impedance matrix Z that translates from the wave amplitudes of the Floquet states , containing the TE/TM amplitudes of the downward ( + , with group velocity in positive z direction) and upward (-) waves f (±) TE∕TM , to the parallel field components in the local coordinate frame. 2 It is thus given by with the q TE and q TM the respective wave number of the fundamental Floquet mode within the air region q (0) 1 , see (4). The corresponding propagation matrix for the Floquet states is with p (±) TE∕TM ∶= exp{± q TE∕TM h} . Using the impedance and the propagation matrix defined above, we thus obtain for the layer transfer matrix The rotation matrix between neighboring layers is with the 2 × 2 identity matrix 1 and the Kronecker product ⊗ . The rotation angle is generally = 2 ∕N ( N ∈ ℕ ), with N = 3 for the hexagonal woodpile with full crystallographic symmetry discussed here. The corresponding unit cell transfer matrix is The eigen decomposition with a diagonal matrix Λ = diag( ) , trivially yields the complex Floquet band structure with the complex wave number The corresponding eigenfield to the Floquet wave number is obtained through re-substituting the field components, given by the -th column of Z −1 .V , into (1).

Chiral woodpile scattering matrix
We here connect the Floquet states belonging to the complex band structure (10) to the scattering parameters reflectivity, transmissivity, and absorbance. We assume a slab of chiral woodpile on a substrate with isotropic (but generally frequency dependent) refractive index n s , and a circularly polarized plane-wave excitation at normal incidence from the top (assumed to be vacuum).
We first introduce the impedance matrix Z 0 (n) that translates from the circularly polarized plane-wave amplitudes within a background with refractive index n to the fields at an the interface at z = 0 in the native coordinate frame, such that The impedance matrix is thus given by Visualization of (q) in (3) for the same parameters as in Fig. 2. The heatmap shows the logarithmic absolute value from zero (dark blue) to large values (red). inherits two branch points (white dots), with corresponding branch cuts (white lines) from the k root functions. The approximate solutions (4) and (5) are indicated by white crosses, while the exact roots are where the ℜ{ } = 0 (solid red line) and the ℑ{ } = 0 (dashed red line) contours cross. Numerically, the first two air roots were found at the green points using a contour integration along the dashed line Naively, the transfer matrix through a finite slab of M unit cells of chiral woodpile is thus given by Using (8) for the unit cell transfer matrix, this expression cannot compute scattering at a semi-infinite slab without substrate, becomes highly inefficient for thick slabs, and most importantly does not reveal the relation to the complex band structure. Using the eigen decomposition (9), however, the slab transfer matrix becomes and for a semi-infinite slab To extract physical meaning, we need to translate these transfer matrices into the corresponding scattering matrices. For this, we first sort the 4 Floquet solutions in (9) with the permutation Π , such that the two waves propagating energy in positive z direction are stored first, and Λ Π ∶= Π.Λ.Π ⊺ and V Π ∶= V.Π ⊺ . 3 We re-express in the sorted basis and can now subdivide all transfer matrices into 2 × 2 sub-blocks that connect downward ( + ) and upward (−) moving amplitudes, respectively: For the scattering matrix we use the most efficient convention [52] to relate the incoming to the outgoing amplitudes: Note that we have included the 2 × 2 phase matrices that transport the incoming amplitudes from the other end of a finite slab to the interface in question for numerical stability [53,54]. Since we are interested in intensities only, these are the identity matrices in the sub-and superstrate, and the respective sub-blocks of Λ M Π within the chiral woodpile.
We can thus express the scattering matrix in terms of the transfer matrix as For the semi-infinite slab, we substitute T inf and P 1 = 1 into (17a) to obtain the reflectivity matrix in the circular polarization basis R = ‖S 11 ‖ 2 . 4 Energy conservation yields the absorptivity A ( ) inf = 1− ∑ � R � for left ( = − ) and right ( = + ) circularly polarized incoming light (from the point of view of the receiver).
For thick finite slabs, the transfer matrix T slab in (12) becomes numerically ill-conditioned. This problem is wellknown [53], and even exists in the fundamental Bragg order here due to the possible strong evanescence of the fundamental Floquet modes (see Sect. 3). As a consequence, the application of (17) to the slab transfer matrix in (12) fails in practice. Instead, the scattering matrix S t at the top interface between the vacuum and the woodpile is well behaved and obtained by substituting into (17). Similarly, the scattering matrix S b between the woodpile and the substrate is obtained by substituting T = Z −1 0 (n s ).V Π , P 1 = Λ M Π,+ , and P 2 = 1 into (17). The scattering matrix through an arbitrary finite slab of chiral woodpile is hence obtained in a numerically well-behaved way through the application of the Redheffer star product [54,57] defined for 2N × 2N matrices C = A ⊛ B with N × N subblocks as The finite slab reflectivity matrix is then R = ‖(S slab ) 11 ‖ 2 , while the transmissivity matrix is T = n s ‖(S slab ) 21 ‖ 2 . As for the semi-infinite woodpile, energy conservation yields (17a) the absorptivity. If the substrate is lossless, we obtain for left ( = − ) and right ( = + ) circularly polarized incoming light. For a lossy substrate, We generally define the (spectral) circular dichroism, both for the slab and the semiinfinite chiral woodpile, as and the spectrally averaged circular dichroism, averaged over the frequency range Ω (from approximately 75THz to 100THz) as

full-wave simulations and materials
The semi-analytical results have been compared to full-wave Maxwell simulations using COMSOL Multiphysics. The frequency domain simulations were performed on a lateral hexagonal unit cell, as shown in Fig. 1a, employing periodic boundary conditions and a tetrahedral finite element mesh with second order tetrahedral elements and a maximum edge length of a/10. We thus simulated a finite woodpile structure for optimized geometrical parameters and 15 unit cells slab thickness. An air domain and a substrate domain of 2c height were added above and below the woodpile, respectively. These vacuum and substrate domains were terminated by periodic ports, in which the system was excited through the vacuum port in the circularly polarized basis. The S-parameters were extracted from the Rayleigh components of the fields at these ports.
The material parameters for the different metals were taken from the refractive index database [58]. We specifically investigated the plasmonic metals summarized in Table 1. As the substrate, we used fused silica glass with n s ≈ 1.5 , silicon with n s ≈ 3.5 , which were approximated as non-dispersive and lossless materials, and the specific metal in question with strongly dispersive refractive index.

Complex band structure of the chiral plasmonic woodpile
To understand the guiding principle behind light propagation in the chiral plasmonic woodpile structures, we start analyzing the evanescent Floquet states calculated by the algorithm introduced in Sect. 2.2. For a broad range of geometrical parameters, which will be discussed in more detail in Sect. 3.2, we find a polarization bandgap in the frequency range of interest Ω between 3 and 4 μ m wavelength: On the one hand, the two fundamental Floquet modes either predominantly couple to left (LCP-philic) or right (RCPphilic) circular polarization. On the other hand, the LCPphilic mode exhibits a band structure with weak dispersion and a small imaginary part, while the RCP-philic mode gives rise to a large bandgap with a relatively large imaginary part and a real part close to the Brillouin zone boundary at ≈ ∕c , as shown in Fig. 4. While the polarization bandgap resembles that of a chiral high index photonic crystal, there are two distinct differences: First, the PC bandgap found here is above the frequency k 0 = ∕c , while in a photonic crystal, it is below that frequency, simplifying top-down fabrication as the corresponding structures for a mid-IR target frequency will be bigger. Second, the RCP-philic mode is not exactly pinned to the Brillouin zone boundary at = ∕c , corresponding to a non-vanishing energy propagation and thus the expected finite energy loss of the strongly evanescent mode within the structure. Similarly, the LCP-philic mode has a small, but non-vanishing imaginary part of , corresponding to a small evanescence or Beer-Lambert-like energy dissipation while propagating through the woodpile PC. As we will demonstrate in Sect. 3.3, the attenuation is underestimated by our approximate theory, as we do not consider the higher-order Floquet modes that concentrate their energy in the metal domain and thus give rise to an additional loss in the PC.
Regarding the circular polarization discrimination, the field polarization of the two Floquet states resides close to the north and south poles of the Poincaré sphere, Table 1 Investigated metals for the chiral plasmonic woodpile, including their permittivity data sources. For each metal, we calculated the maximum of the average CD over the spectral range Ω spanning from 75 to 100 THz, obtained for a semi-infinite PC of fixed layer pitch d = 2.4 μm , and variable volume fill fraction of the metal and height of a single layer h (cf. Fig. 1). The average CD at the optimized position decreases for a finite slab of 15 unit cells on a substrate made of glass ( n s = 1.5 ), silicon ( n s = 3.5 ), or the metal in question (m) respectively. We quantify the circular dichroism of the Floquet states with the circular dichroism index [12,59] that considers the relative difference in coupling to the two circularly polarized plane waves and the Floquet field and ranges from 0 (no difference) to ±1 (maximum difference). The total in-coupling from the vacuum is similarly quantified through the coupling index that ranges from 0 (no in-coupling) to 1 (maximum incoupling). The individual RCP and LCP couplings in these expressions are approximated by where V (i) Π ( i = 1, 2 ) is the ith sorted eigenvector of the unit cell transfer matrix, that is the ith column of V Π , which contains the lateral electromagnetic fields of the corresponding Floquet state on the C 6 2 axis at the top of the unit cell. As this definition considers both the electric and magnetic fields, it takes the impedance match into account.
We find that the circular dichroism index for the platinum structure corresponding to the band structure in Fig. 4 at wavelengths 0 = (3, 3.5, 4) μm is C = (−0.84, −0.83, −0.81) for the LCP-philic mode and C = (0.98, 0.94, 0.89) for the RCP-philic mode. In other words, incoming LCP light couples almost exclusively into the blue, weakly attenuated propagating Floquet state, while RCP light couples predominantly into the red, strongly evanescent state in Fig. 4a, justifying the color assignment. Further, the coupling index is = (0.99, 0.99, 0.98) into the LCP-philic mode, and = (0.67, 0.52, 0.52) into the RCP-philic mode for 0 = (3, 3.5, 4) μm . Even though the approximate theory is expected to overestimate the coupling, this suggests that almost all of the incoming LCP light is transmitted into the woodpile PC and absorbed through Beer attenuation within. On the other hand, a considerable amount of RCP light is expected to couple into the red, strongly evanescent Floquet mode. Since this mode, however, propagates almost no energy into the structure within the bandgap region, RCP light is mainly reflected with little absorption.

Optimizing broadband circular dichroism
While the qualitative chiro-optical behavior of the chiral woodpile PC is well understood within the band structure picture discussed in Sect. 3.1, we here perform a more quantitative analysis calculating the scattering matrix and the circular dichroism as outlined in Sect. 2.3. Next to revealing the main physical scattering mechanism, the semi-analytical Floquet mode algorithm comes at a very low numerical cost. The scattering parameters for one wavelength are calculated within 2 ms (in a non-optimized Python implementation), compared to 5 min for the full-wave simulations (see Sect. 2.4 for details). This efficiency enables a brute-force scan over geometrical parameters for different materials to extract the spectrally averaged circular dichroism ⟨CD⟩ , defined in (20). We find, that a strong CD can be generally observed for lattice constants below, but on the order of 3 μm , and have therefore fixed the single layer pitch to d = 2.4 μm ( a ≈ 2.77 μm ) for the remainder of this manuscript.
For the four metals listed in Table 1, we have calculated the average CD of the semi-infinite metal-air PC with varying thickness of the in-plane air region d 1 ∈ [1, 2] μm and the layer height h ∈ [0.3, 1.5] μm . Note that the metal fill fraction and the chiral pitch are hence varied according to = 1−d 1 ∕d and c = 3h . The results illustrated in Fig. 5 reveal a broad region within the parameter space, where a substantial CD can be observed, largely independent of the chosen metal. This region between 600 and 800 nm layer heights, where the spectrally averaged CD is above 0.5, spans almost all d 1 , with a slight increase towards larger thicknesses. While aluminum and tungsten yield the strongest CD for the semi-infinite slab (cf. Table 1), their more strongly perfect plasmonic nature ( ℑ{n} ≫ ℜ{n} ) gives rise to very little field penetration into the metal region. This in turn makes them less efficient in more realistic finite slabs, since the attenuation of the LCP-philic mode (Sect. 3.1) and thus LCP absorption is much weaker. As discussed in Sect. 3.3, attenuation is, however, underestimated by our theoretical approximation. The chiral PC is expected to yield a substantial spectrally averaged CD, (20), of 70% for all metals, insensitive to fabrication imperfections. A specific metal alongside geometrical parameters within a broad region may hence be chosen depending on the limitations of a specific fabrication routine.
In summary, the CD mainly depends on d 1 and h (and not explicitly on a) and is strongest within a specific region in Fig. 5, which is independent of the metal in question. To better understand this behavior, we approximate the chiral woodpile even further by replacing it by a plasmonic version of a semi-continuous Bouligand structure, as found in liquid crystals and biological systems [60]. Since we expect the electromagnetic field to rotate with the smaller 60 • rotation of the C 6 screw rotation, and not the 120 • of the C 3 rotation used in Sect. 2.2 and illustrated in Fig. 1, we build the Bouligand structure such that it rotates in the opposite direction and has a pitch of 2c. The anisotropic lateral (x-y) permittivity matrix depending on TE∕TM ∶= q TE∕TM ∕k 0 at height z is in the Bouligand picture expressed by: The monochromatic Maxwell equations can be solved with the ansatz for the lateral electric field (with E z = 0 ). With the approximate solutions for the Floquet wave numbers (4), this procedure yields the following 2D quadratic eigenproblem in : with G 1 ∶= ∕d 1 . One can immediately see, that this approximate Bouligand equation only depends on the chiral pitch c (or layer height h) and thickness of the air region d 1 . and G ∶= ∕c. Evidently, − is purely imaginary, as expected in the center of the bandgap. It connects the R solutions and is linearly polarized with polarization direction between x and y. On the other hand, + belongs to the L branch and is left elliptically polarized, quickly approaching circular polarization if h ≪ 2 3 d 1 . The Bouligand model thus identifies three main contributors that limit the region of high average spectral CD. First, the R bandgap needs to reside within the spectral range of interest. To keep the lower band edge outside the spectral region of interest, we require G < 2 ∕ l , with l = 4 μm in our case. This yields h > l ∕6 ≈ 0.67 μm , shown as brown line in Fig. 5, which explains the horizontal border at the bottom of the high CD domain. Similarly, keeping the upper bandgap edge outside the spectral region of interest requires (where the radicand is positive), shown as a red line in Fig. 5. Above this line, the red region in the heatmap indicates a sign change in the average CD. This change is caused by a dichroic color switch [61], meaning that there is an additional bandgap of opposite optical chirality at higher frequencies, which also exists in chiral woodpile photonic crystals [38] and moves into the spectral region of interest. Finally, the polarization of the L mode needs to resemble LCP polarization, i.e. it needs to be as close as possible to the pole of the Poincaré sphere, to yield a strong CD in the bandgap. This implies that we need to be as far as possible under the white h = 2 3 d 1 line in Fig. 5, explaining the less sharp positively sloped upper termination of the high CD region. While the Bouligand model can thus explain the general shape of the high CD region, it predicts the red bandgap mode to be linearly polarized in contrast to the hexagonal chiral woodpile bandgap modes, which are clearly right circularly polarized as demonstrated by the CD index in Sec. 3.1.

Comparison to full-wave simulations
Considering the crude approximations made, the theory predicts the full spectral CD of all metals well. To demonstrate this, we have calculated the RCP and LCP absorption through full-wave simulations as described in Sec. 2.4. The results for an iron-air woodpile PC on a glass substrate with d = 2.4 μm , = 0.2 , and h 1 = 0.7 μm , close to the optimal parameters of the semi-infinite slab, are shown in Fig. 6. Clearly, the position of the bandgap is predicted well, although the full-wave results seem to be slightly blueshifted compared to the theory. Generally, the theory underestimates the absorption for both polarizations within the bandgap region. This behavior is expected, as we do not consider the single layer higher order fields that reside mainly in the metal domain. These lead to stronger absorption both in the evanescent Floquet mode (red spectrum), and the propagating mode (blue spectrum). An additional indication that absorption and hence attenuation in the propagating Floquet mode is increased in the simulations is the absence of the Fabry-Pérot interference pattern that is clearly visible in the theoretical spectrum. At high wavelengths above 4 μm on the other hand, where iron begins to act more and more like a perfect electrical conductor, the inclusion of the higher order metallic modes is expected to lead to an underestimated impedance mismatch between the incoming vacuum field and the Floquet states within. The theory, therefore, underestimates the reflection and overestimates the absorption for both polarizations alike.
While the theory cannot accurately predict the spectra, we have thus demonstrated that all qualitative physical predictions made in the preceding chapters are accurate and can be used to tailor a chiral woodpile PC to a specific application and fabrication procedure. These general principles can also help to further fine-tune the geometry to improve the nonapproximated CD. For example, the position of the bandgap in Fig. 6 seems to be too far to the left. This suggests increasing the layer height h slightly, which is predicted to red-shift the right band edge towards 4 μ m, while only weakly affecting the left edge at 3 μ m. We demonstrate this effect by simulating the same structure with an increased h 2 = 0.72 μm , shown as dotted line in Fig. 6.

Conclusion
In conclusion, we have identified the underlying physical principles of broadband circular dichroism in chiral plasmonic woodpile structures employing an approximate evanescent Floquet mode picture. Focusing on the transparency window of the atmosphere between 3 and 4 μ m wavelength in the mid-infrared frequency region, we have found a broad circular polarization bandgap that exists within a large region of geometrical parameters and for a number of different metals used. Employing a semi-continuous Bouligand model, we extracted general design principles that predict the approximate size and shape of the region in the geometrical parameter space, where large broadband circular dichroism is expected.
On the one hand, our findings demonstrate that evanescent Floquet modes and the associated complex bandstructure form an invaluable tool to understand scattering at and wave propagation within slab-like plasmonic crystals, which combine interference-dominated physics known from photonic crystals and material-dispersion induced effects known from classical metamaterials. On the other hand, we provide a pathway to design a broadband, highly efficient circularly polarized thermal source in the mid-infrared region by application of Kirchhoff's law. A large broadband CD can be engineered within a predictable, massive region in the geometrical parameter space. This region encompasses a variable aspect ratio of the metal bars, which ranges from a ratio of approximately 1 : 2 to 2 : 1. The typical lattice constant is smaller, but comparable to the wavelength of the light, approximately twice as big as in high-index dielectric structures, making top-down fabrication more feasible for the envisioned mid-infrared window between 3 and 4 μ m wavelength.
Indeed, a number of fabrication routines have been reported to yield woodpile structures for geometrical parameters within the predicted high-CD region, ranging from layer-by-layer manufacturing [62,63] to two-photon lithography techniques [64]. Recent advances in two-photon lithography make it possible to directly produce metallic structures [65], while all other methods can produce an inverse mold on a conducting substrate (for example indium tin oxide). In both cases, an electro-deposition routine [66] can replicate the woodpile PC in a metal of choice. Standard direct laser writing on the other hand typically produces a polymeric woodpile structure, for which an electroless plating routine can be employed [67,68].

Conflict of interest Matthias Saba is an Editorial Board Member of Applied Physics A. All other authors declare no conflict of interest.
Open Access This article is licensed under a Creative Commons Attribution 4.0 International License, which permits use, sharing, adaptation, distribution and reproduction in any medium or format, as long as you give appropriate credit to the original author(s) and the source, provide a link to the Creative Commons licence, and indicate if changes were made. The images or other third party material in this article are Fig. 6 Theoretically predicted (solid lines) and simulated (dashed lines) absorption spectra under LCP (blue) and RCP (red) illumination for an iron-air chiral woodpile PC on a glass substrate ( n s = 1.5 ) with d = 2.4 μm , = 0.2 , and h 1 = 0.7 μm . The right band edge red-shifts in the simulated absorption spectra if the layer height is slightly increased to h 2 = 0.72 μm (dotted lines) included in the article's Creative Commons licence, unless indicated otherwise in a credit line to the material. If material is not included in the article's Creative Commons licence and your intended use is not permitted by statutory regulation or exceeds the permitted use, you will need to obtain permission directly from the copyright holder. To view a copy of this licence, visit http:// creat iveco mmons. org/ licen ses/ by/4. 0/.