Alignment-insensitive bilayer THz metasurface absorbers exceeding 100% bandwidth

: Metamaterial absorbers have been a topic of considerable interest in recent years, with a particular focus on Terahertz (THz) frequencies due to many natural materials having a weak interaction with THz light. Great eﬀorts have aimed to expand such THz absorbers to cover a wide bandwidth whilst also being highly eﬃcient. However, many of these require cascaded or stacked multilayer resonant elements, where even a small deviation in the alignment between layers is extremely detrimental to the performance. Here, we propose a bilayer metasurface absorber (thickness ∼ λ /6) that is immune to such layer misalignments capable of exceeding a fractional bandwidth (FWHM) of 100% of the central frequency. The design works due to a novel absorption mechanism based on Salisbury Screen and anti-reﬂection absorption mechanisms, using fractal cross absorbers to expand the bandwidth. Our work is of particular beneﬁt to developing devices which require ultra-wide bandwidth, such as bolometric sensing and planar blackbody absorbers, with the extremely robust absorption responses being unaﬀected by any misalignments between layers – a limiting factor of previous absorbers. work

A particular emphasis has been placed on advancing these devices in the THz regime, which is classed as being the frequency range of 0.1-10 THz (wavelength 30 µm -3 mm), due to the low interactivity of THz radiation with normal materials; many applications can be utilized from THz radiation, such as spectroscopy of several interesting materials (including illicit drugs and explosives) which have characteristic "fingerprints", imaging of biological tissues due to its ability to penetrate visibly opaque media (and THz radiation being non-ionizing and biologically safe), and other non-destructive and security related examples [32][33][34][35][36][37][38][39][40][41][42][43][44]. One device which would be of particularly important value is the THz blackbody perfect absorber (emitter), which would allow complete absorption (emission) of light in the whole THz spectrum. Extensive works have been published on the development of THz absorbers, with the aim of extending the bandwidth as large as possible whilst also being highly efficient -the key requirements for such a perfect blackbody absorber. However, many of these absorbers suffer from certain setbacks -the primary methods of achieving a wide absorption bandwidth involves either cascading resonant elements in the vertical or horizontal directions. The former approach has yielded impressive results [45], but suffers from fabrication issues -small deviations or imperfections in the alignment of the elements results in a much poorer performance, as was investigated by Chen et al [46], where it was found that a λ/10 misalignment reduces the absorbance from near 100% to 30%. The latter approach of horizontally arranging resonant elements allows much more freedom and simplicity, by typically allowing a single-layer and fabrication step, but does not yet match the efficiencies and bandwidths achieved by stacked absorbers [47,48]. This is due to the fact that for different resonant elements they typically require different dielectric spacer thicknesses to optimize the absorption through ground-plane coupling.
Recently, we investigated horizontally arranged fractal elements as a means of achieving wide (octave-spanning) absorption bandwidth [48]. This result was based on a novel absorption mechanism involving Salisbury Screen (SS) and anti-reflection (AR) responses, allowing the absorbance to also be simultaneously large. Many previous approaches of this kind used much thinner dielectric spacers, and therefore relied on the typical inductive coupling of the elements and ground-plane to produce the absorption [17], which is found to be detrimental for combining multiple resonant elements together due to the destructive interference of the induced currents. Instead, we showed that a slightly thicker spacer yields a much more broadband absorption in both experiment and simulation, as there are negligible ground-plane currents present, surpassing many other planar absorbers in performance.
In this work, we utilize the combined SS/AR mechanisms of resonant fractal cross structures alongside the introduction of an embedded fractal cross layer within the dielectric spacer ( Fig. 1(c)) to achieve excellent absorption bandwidth. We exhibit a record fractional bandwidth (FWHM/center frequency) of 103.7% (equivalent to almost 5 THz bandwidth) and an average absorption efficiency of 87.4%, with the minima in the absorption spectra never dipping below 80%. This bilayer design has the additional benefit of being insensitive to layer misalignments, where a complete misalignment of the layers results in only a 9% variation in the fractional bandwidth and a 3% decrease in the average absorption. The design boasts a 53% improvement in the bandwidth over our previous single layer design (which had a fractional bandwidth of 68%) [48], therefore exhibiting the broadest bandwidth planar THz absorber to date. Additionally, the design is well-performing even with incident angles up to 45°, where the bandwidth is unaffected and the absorbance only slightly decreases to ∼60% at lower frequencies. Such a design will be of particular benefit to the development of low form-factor blackbody THz absorbers and bolometric sensors, where the bilayer alignment insensitivity will be particularly attractive for commercialization.

Characteristics of fractal crosses and bilayer design
The design used in this study relies on the well-known metal-insulator-metal (MIM) absorber design, as shown in Fig. 1(a), whereby resonant elements (termed fractal, length d) are placed above a metallic ground-plane (thickness t) separated by a dielectric spacer (in this case polyimide, thickness h). The design allows the incident THz waves to undergo a Fabry-Pérot type resonance effect, where a minimization of the reflected light R(ω) corresponds to a maximum in absorption via the equation A(ω) = 1 -T(ω) -R(ω), where T(ω) = 0 due to the presence of the metal ground-plane. Tantamount to our design is the use of resonant crosses undergoing fractal order increases, as shown in Fig. 1(b). Increasing a shape's fractal order relies on the addition of a self-similar resonator shape onto itself, which works to increase the effective electric length of the original resonator to produce a red-shifted resonance frequency [49,50]. Therefore in Fig. 1(b) we attribute the increase in fractal order from fractal cross 1 to 4 with lowering of the resonance frequency due to the increasing fill factor. In our device, a factor of 1/3 is used to scale the first cross onto itself, to become the fractal 2 cross, and then by a factor of 1/2 for the following crosses. The values given are then d1 = 20 µm, d2 = 6.67 µm, d3 = 3.33 µm, and d4 = 1.67 µm, where the periodicity of the crosses p = 32 µm. The thickness of all crosses are 0.1 µm (t/2) and the ground-plane thickness was set as 0.2 µm (t), sufficiently thick so as to prevent light transmission. The width of all resonator 'wires' are set as 0.5 µm, as any larger would potentially lead to touching arms for fractal level 4. Both the crosses and ground-plane are chosen as gold, with a conductivity in the THz regime as σ = 4 × 10 7 S/m. Our choice of dielectric spacer is polyimide (PI2545) with a refractive index (obtained by best matching to experimental results in [48]) of n = 1.71 + 0.08i, and a thickness h = 11 µm.
For our bilayer design, smaller crosses were used for the bottom layer embedded within the polyimide spacer, shown in Fig. 1(c). The size difference is due to the fact that the surrounding environment of these crosses is completely polyimide as opposed to the surface crosses which have polyimide substrate and air superstrate -the increase in the refractive index of the dielectric environment then results in a shift in the resonance frequency of the fractal crosses. Combined with the fact that these embedded crosses are positioned at a different height above the groundplane compared to the top layer, it was critical to reduce their sizes in order to best optimize absorption at higher frequencies. The size of the bottom layer crosses are therefore scaled down by a factor of ∼ 1/ √ 2. For the best performance, the embedded crosses were positioned a distance l = 4 µm below the top layer (polyimide/air interface), corresponding to a distance of 7 µm above the ground-plane. The supercell designs differ between the top and bottom (embedded) fractal cross layers, as is shown in Fig. 1(d). The maximum fractal order which is used in the bottom layer only goes up to a level of 3, compared to the top layer which goes up to fractal level 4. We do this for two reasons, with the first being that because of the limited width of the original (top layer) fractal cross wires already being relatively small (0.5 µm), the resulting width of the scaled down crosses would be even further reduced, and be both challenging to fabricate and to simulate (decreased mesh size and increased computation time) if a fractal level of 4 was utilized. The second reason relates to the fractal level 1 bottom cross, which is the highest frequency resonator used in the design; to ascertain a good absorption efficiency in the broadband device at this highest frequency, it was necessary to use two fractal level 1 crosses in the bottom layer (shown in Fig. 1(d), right) compared to only one each of the fractal level 2 and 3 crosses. Comparatively, when using only one fractal level 1 cross for the bottom layer in preliminary designs, the absorption at the higher frequencies was significantly reduced which we attribute to the low cross section and fill factor.
In order to assess the design functionality of our bilayer absorber, we performed numerical simulations using the commercially available CST Microwave Studio software. All simulations were carried out in the Frequency Domain solver, where periodic boundary conditions were applied in the x and y directions (in-plane, parallel to absorber) with the boundary conditions in the z direction (out-of-plane, perpendicular to absorber surface) being open for the port signals. Using the dimensions given in Fig. 1, arrays of individual fractal crosses were simulated and their corresponding resonance frequencies obtained, as given in Fig. 2. Figure 2(a) shows the individual cross spectra for the top layer, whilst Fig. 2(b) shows the spectra for the bottom layer crosses. We see in Fig. 2(a) that a distinct red-shift of the resonance frequency occurs with increasing fractal order, with maximum absorption peaks located at 5.3, 3.33, 2.62 and 2.47 THz for fractal level 1-4 respectively. It is noticeable that the fractal level 2 cross (red curve) has 2 resonant peaks situated at 3.33 THz and 4.62 THz. We attribute the resonant peak at 4.62 THz to be related to the anti-reflection response of the polyimide itself, as opposed to the fractal cross structure. This is confirmed by removing all fractal structures in the simulation and obtaining the dashed orange curve (labelled Polyimide Only) with a clear peak occurring at 4.65 THzvery close to the 4.62 THz of the fractal level 2 cross. It is also noticeable from Fig. 2(a) that the other fractal structures have secondary peaks (shoulders) in their spectra also located around the (∼4.65 THz) polyimide anti-reflection resonance, thereby providing further confirmation. Arrays of scaled-down individual crosses embedded within the polyimide were investigated and their resonance spectra shown in Fig. 2(b). The results show similar trends to the larger crosses (as in Fig. 2(a)), albeit occurring at higher frequencies. The peak resonance frequencies occur at 5.50, 4.32, and 3.49 THz for fractal level crosses 1-3 respectively. We also show the dashed orange anti-reflection curve for reference. Due to these crosses being embedded within the polyimide and not at the boundary, the anti-reflection response has less of a pronounced effect on the fractal cross spectra, evident by the lack of secondary peaks or shoulders. From Figs. 2(a) and 2(b), we can therefore deduce in principle that the combined bandwidth of these absorbers span from 2.47 THz to 5.50 THz (at the peak wavelength).
Following on from this, we combined the corresponding layer crosses into the supercell arrangements as given in Fig. 1(d) and examined their respective spectra. The results of these are given in Fig. 3. We see that for the combined top fractal crosses (green curve) there is a relatively consistent high absorption exceeding 70% from 2.58 to 5.07 THz, almost in agreement with the frequencies spanned for the individual cross-array resonance frequencies (2.47 to 5.3 THz). We attribute the variation in both bandwidth and absorption magnitude (at higher frequencies) to cross-cross coupling affecting the local resonant frequencies of the crosses, resulting in a difference compared to when the crosses are independently arranged. Similarly, the combined bottom layer fractal crosses (magenta curve) also show a consistently high broadband absorption exceeding 70% from 3.48 to 5.27 THz. Again, this frequency range almost agrees with the resonance frequency span of the individual crosses (3.49 to 5.50 THz) in the bottom layer, with a slight deviation again due to cross-cross coupling. However, because the bottom layer crosses are reduced in size, yet maintain the same cross-cross (center-center) distance (periodicity), coupling is less pronounced compared to the top layer arrangement. Our final design has both the top and bottom fractal layers combined into one bilayer absorber, as previously shown in Fig. 1(c). The resulting spectra of the bilayer design is given as the blue curve in Fig. 3. It is clear that this combined spectra almost perfectly aligns with the spectral responses of the individual top and bottom layers (green and magenta curves) in Fig. 3, thereby giving a very broadband and highly-absorbing response. The spectra shows a very flat absorption response over 90% spanning from (peaks) 2.63 to 5.12 THz. The highest absorption value which occurs is 99.6%, resulting in a "half-max" of 49.8%. Using this, we calculate the Full-Width at Half-Max (FWHM) ∆f of our bilayer absorber to be 4.22 THz (between the frequencies of 2.35 and 6.57 THz). Correspondingly, we can calculate the central frequency f c as 4.46 THz.
Using both ∆f and f c , we can then deduce the fractional bandwidth of our absorber [48] using ∆f / f c , which gives us a value of 94.7%. As a direct comparison, our previous single layer fractal design [48] which was the broadest bandwidth THz absorber at the time of publication achieved a fractional bandwidth value of 68% (corresponding to a bandwidth greater than an optical octave). The proposed bilayer design here improves upon this with almost a 40% increase in bandwidth. One limitation of proposed bilayer designs in the past, however, has been ensuring precise alignment (λ/10) between layers to obtain the optimum absorption characteristics.

Effect of misalignment between layers
In this section we examine the aforementioned effect of misaligning the top and bottom layers laterally. This is of considerable importance for fabrication of such absorber devices, where fabrication tolerances and alignment may not be exactly as designed with slight deviations potentially leading to performance issues. As was mentioned in Section 1, work carried out by Chen et al. [46] showed that a very small deviation (λ/10) in the cross positions on different layers resulted in a decrease of absorbance by 70%. Here, we now show that using the modified broadband absorber design we obtain a device which is not only insensitive to absorption decrease with layer misalignment, but actually performs slightly better when the crosses are completely misaligned. Figure 4 gives the spectra for both x and y polarized light normally incident onto the bilayer absorber, for four designs with different alignmentscomplete alignment ( Fig. 4(a)), x misalignment ( Fig. 4(b)), y misalignment (Fig. 4(c)), and complete misalignment in both x and y directions ( Fig. 4(d)). On first glance it would be assumed that both designs, Fig. 4(b) and 4(c), should give identical responses, seeing as both x and y polarizations are used, but this is not the case. The difference is due to the fact that the top layer appears chiral resulting in differing responses to x and y polarizations. The bottom layer is, however, both mirror and rotationally symmetric (C4), and will provide identical responses to orthogonal polarizations. The resultant spectra of the two layers will therefore behave differently for both x and y polarizations/misalignments.
It is observable that all spectra in Figs. 4(a)-4(d) have absorption efficiencies exceeding 80% over large bandwidths, irrespective of either polarization or layer alignment. Even more interestingly is the bandwidth improvement when going from completely aligned layers (Fig. 4(a), blue curves) to completely misaligned layers (Fig. 4(d), black curves), which is in contradiction with the previous literature [46]. To better assess the differences exhibited by the metamaterial absorbers with differing alignments, we quantitatively compare all four designs below in Table 1. Notable figures of merit included in the table are both min and max peak absorption values, FWHM, fractional bandwidth, and the average absorption. We see that for all of the absorbers, the lowest peak absorption (the lowest absorption value which occurs in the absorption plateau above 80%) is 82.4%, whilst the highest peak absorption reaches at least 99.6% for all absorbers. As seen in Figs. 4(c) and 4(d), there appears to be a small sharp dip situated at ∼4.76 THz. We believe this to be a diffraction-type effect, owing to the periodicity of the scatterers. Typical diffraction would occur at a wavelength corresponding to the periodicity of the scatterers, which in this case is 32 µm -this is equivalent to a frequency of 9.4 THz which is far above our frequency range of interest, and so can be ruled out. However, because the supercell has a periodicity of 64 µm (which is the distance between identical crosses of the same resonant frequency), there could exist some imperfect diffraction-type effect where diffraction can occur at 64 µm wavelength (equivalent to 4.69 THz, close to the dip at 4.76 THz) -but due to the presence of de-tuned resonators (namely crosses with a broadband resonance which is offset in frequency) positioned at 32 µm in between, typical diffraction is prevented. Another phenomenon to mention is the Wood's Anomaly (WA) effect [48,51,52], which occurs at a wavelength corresponding to λ WA = n · P, where n is the refractive index of the substrate and P is the distance between resonators. In this case, n ≈ 1.71 and P = 32 µm, which gives a wavelength of λ WA = 5.47 µm equivalent to 5.48 THz, also very close to the spectral dip of 4.76 THz in Figs. 4(c) and 4(d). We feel that these diffraction-type effects have an almost negligible influence on our result though, and so we do not investigate this further.  It can sometimes be misleading when the literature states the FWHM as the key figure of merit (FOM) defining the bandwidth performance, as a larger FWHM may lead the reader to assume a device is more broadband simply by being at a region of higher frequency. Instead, we propose an alternative figure of merit for assessing how broadband an absorber is, by calculating the fractional bandwidth, which was stated in section 2 as being the ratio between the FWHM and central frequency of the band. It is a more reliable FOM than the FWHM, as it gives a bandwidth percentage normalized to the frequency range of interest, regardless of frequency, and is especially useful when comparing two different absorber designs at different frequency regimes. From this, we can directly compare the four different absorber designs shown in Fig. 4 and state the fractional bandwidths in Table 1. It is observable that the lowest absorption bandwidth occurs for y polarized light incident on the aligned design ( Fig. 4(a)) with a fractional bandwidth of 94.7%, whilst the best performing absorber achieves a value of 103.7% when the layers are totally misaligned (Fig. 4(d)). To the best of our knowledge, this is the best performing THz bandwidth achieved to date for a bilayer metasurface absorber (superseded only by more complex multilayer Fabry-Pérot type absorber designs [53,54]), and achieves a 53% improvement upon our previous single layer absorber design (which achieved an experimental fractional bandwidth of 68% [48]). It may be noticeable that the lowest and highest fractional bandwidths listed in the table almost coincide with the lowest and highest FWHM, of 4.224 THz and 4.822 THz respectively. The highest FWHM actually occurs for the y-misaligned design, with a FWHM of 4.845 THz -slightly higher than the totally misaligned design -but only has a fractional bandwidth of 102.4%, further reinforcing the necessity of using the fractional bandwidth instead of the FWHM as a figure of merit for bandwidth performance. If we assume unpolarized light and do not choose to align the top and bottom layers, we could hope to expect an average fractional bandwidth of 100.6% when taking the average of all values given in the table. This is a very impressive result, whereby we can conclude that the average performance of our device is an absorber exceeding 100% bandwidth regardless of the alignment or polarization. Such a device would be of particular benefit to commercialization as it would not involve stringent alignment tolerances and therefore much more relaxed fabrication processes could be used.
The absorption magnitude is an equally significant FOM for a broadband absorber. All absorber designs shown have a minimum peak absorbance of at least 82.4% (as stated previously), but this value is only defined for a single specific frequency. A more appropriate FOM would be the average absorption of the device -which we define as being the average absorption value between the minimum and maximum frequencies of the FWHM. As the peak absorption value for all designs is essentially 100% (lowest being 99.6%), the half-maximum value of the FWHM is 50%, and therefore the average absorption values, as given in Table 1, are the averaged sum of all frequency dependent absorption data above 50%. From this, we can see that the average absorption above the FWHM is at least 83.5% (for the y misaligned design), and reaches 87.4% for the totally misaligned design. Therefore, we can state that the totally misaligned design has both the broadest bandwidth (103.7%) and highest average absorption (87.4%) of the four designs, which is significant as it opposes the typical approach that multilayer metasurface absorbers must be aligned to achieve the best result.

Electric-field profiles of top and bottom layers
In order to better understand the phenomena occurring within our absorber design, we performed numerical simulations to illustrate the electric-field profiles of the two layers at a range of frequencies. These frequencies of interest are those which correspond to the peaks in the x polarized spectrum for the totally-misaligned design in Fig. 4(d). Numerical simulations for the electric-field profiles were performed in Lumerical FDTD software, and so resulted in a slight change in the peak positions of the broadband spectra for x polarization incident on the totally misaligned design. However, the variations from the CST simulation results are minimal and therefore we assume that the underlying physical phenomena taking place are still valid. A total of six frequencies were chosen to perform the simulations for obtaining the 2D electric-field profiles given in Fig. 5, below. The six frequencies chosen were 2.72, 3.45, 4.64, 5.08, 5.53, and 6.59 THz which correspond to Peaks I -VI (of the absorption spectra for the misaligned configuration), respectively. Figs. 5(a)-5(f) correspond to the 2D E-fields of the top absorbing layer for these six frequencies, whilst Figs. 5(g)-5(l) correspond to the 2D E-fields of the bottom absorbing layer crosses for the six frequencies -corresponding frequencies are shown above the sub-figures in Fig. 5.
A general trend which can be seen for both top and bottom layers are the enhanced 2D E-field profiles which occur on the higher level fractal crosses for lower frequencies, which shift to the low level fractal crosses as the frequency increases. The arrangement of the crosses for the top layer is the same as that in Fig. 1(d) (left), namely the fractal level increases from level 1 to level 4 in a clockwise direction from the top-left to the bottom-left, respectively. The bottom layer has the same cross arrangement as Fig. 1(d) (right), where the unit cell is centered on the level 1 cross. Due to the fact that there are two level 1 crosses per unit cell in our design, the level 2 and level 3 crosses are visible on the left/right and top/bottom borders, respectively, whilst the second level 1 cross occupies the four corners of the unit cell. Whilst this bottom layer configuration is not exactly that as given in Fig. 4(d), the spectra are almost indistinguishable due to the polarisation insensitivity of the bottom layer, showing the robustness of our design.
For the lowest of the six chosen frequencies (Peak I = 2.72 THz), the closest cross resonances are those of the top layer level 3 and level 4 crosses which resonate at 2.62 and 2.47 THz, respectively (as seen in Fig. 2(a)). Indeed, we see in the E-field profile of Fig. 5(a) that the level 3 cross has the strongest field profile, which corresponds to its resonant frequency of 2.62 THz being close to 2.72 THz, and explains why the level 4 fractal cross (2.47 THz) in Fig. 5(a) does not have a stronger E-field profile. It can be justified that whilst the individual lowest resonance frequency is 2.47 THz for the isolated fractal level 4 crosses, it does not occur in the broadband spectra where 2.72 THz is the lowest frequency peak. As shown in Fig. 2(a), the resonances of both the fractal level 3 and level 4 crosses are extremely close together, which, when combined together in the supercell arrangement of Fig. 1(d) give a different spectral response as in Fig. 3 (green curve). This is evidenced by the fact that the lowest individual resonance for the fractal level 4 cross occurs at 2.47 THz (Fig. 2(a)), whilst the lowest peak frequency which occurs for the combined top layer crosses (Fig. 3) is increased to 2.58 THz. We attribute this to a hybridization of the very close resonances between the fractal level 3 and level 4 crosses, which may have some coupling interactions due to being in close proximity to each other (12 µm distance between the cross extremities), thereby resulting in a slight blue-shift of the minimum frequency.
Due to the broadband nature of the fractal crosses used in our design, we also see that at 2.72 THz all of the top layer crosses undergo some excitation and exhibit enhanced E-field profiles. As the frequency increases to 3.45 THz in Fig. 5(b) (Peak II), the E-field excitation primarily shifts to the fractal level 1 and level 2 crosses. Again, when the frequency increases to 4.64 THz in Fig. 5(c) (Peak III), the excitation shifts even more to the fractal level 1 cross, as expected. For higher frequencies (Peaks IV-VI in Figs. 5(d)-5(f)) the excitation seems to gradually diminish on all crosses -this can be attributed to the differences in the combined layer bandwidths which are shown in Fig. 3, and therefore these higher frequencies primarily corresponding to the bottom layer (smaller) crosses. This can be seen in Fig. 5(j) (5.08 THz) and Fig. 5(k) (5.56 THz) where the cross excitations are maximized and E-field profiles are strongest, with 5.08 THz being the approximate mid-point of the combined bottom layer spectrum given in Fig. 3 (magenta curve). Although the bottom layer crosses have higher resonant frequencies than the top layer crosses, some frequency overlap still occurs (the green and magenta curves in Fig. 3), which therefore explains why there are still some strong E-field excitations occurring for the bottom layer crosses at the lower frequencies in Figs. 5(h) and 5(i). To verify the absorption mechanism taking place we performed numerical simulations in Lumerical FDTD to view the absolute E-field profiles (|E|) in the y-z planes, in order to see how the energy is localized at resonance. This helps us to envision where absorption primarily occurs. A totally misaligned arrangement was chosen, corresponding to the arrangement in Fig. 5, so that we could distinguish the resonances taking place in different areas of the simulation -these are shown in Fig. 6. We chose y-z cuts at three x positions, x = −16 µm, 0 µm, and + 16 µm, corresponding to cuts through both the top and bottom layer crosses (we choose the center of a supercell as being 0 µm). A schematic is given in Fig. 6(a), where the eye signifies the direction that we are looking at the cut-planes for reference to Figs. 6(b) and 6(c). To get a better insight into the absorption mechanism, we chose two frequencies to examine, 2.72 THz (Fig. 6(b)) and 5.08 THz (Fig. 6(c)), which are close to the "shoulders" of the broadband plateau given in Fig. 4(d). As seen in Fig. 6(b), there are very high localizations of E-fields particularly for the fractal level 4 cross (at x = −16 µm) and fractal level 3 cross (x = +16 µm). From Fig. 2, we know that the resonance of the fractal level 3 cross is 2.62 THz, which is very close to the investigated frequency of 2.72 THz in Fig. 6(b), and confirms the reasoning for the level 3 cross having the highest localization of electric field. It is noticeable that the E-fields are nearly entirely localized around the vicinity of the metallic crosses, and quickly decay towards the ground-plane -we can therefore deduce that absorption primarily takes place near the resonant crosses and the ground-plane does not play a significant role in absorption, as is the case for the previous generation of narrow-band absorbers [45]. Similarly, for the E-field plots at 5.08 THz in Fig. 6(c) we expect the smaller crosses (in the bottom layer) to show more prominent E-field localization due to the higher resonant frequencies (as given in Fig. 2). This appears to be the case, where the bottom-layer fractal level 1 cross shows very high E-field localization, also in agreement with Fig. 2(b) where the level 1 cross has its resonance around 5.50 THz. These results are also very well matched with the x-y E-field profile plots given in Fig. 5, where Figs. 5(a) and 5(g) correspond to Fig. 6(b) at 2.72 THz, and Figs. 5(d) and 5(j) correspond to Fig. 6(c) at 5.08 THz, respectively.

Absorption response to incidence angle
In order to assess the potential robustness of our design, we performed numerical simulations using Lumerical FDTD software to examine the effect of oblique angle and polarization of the light incident onto our absorbers. The design used for these simulations was the full bilayer arrangement, with the bottom layer arranged as for the Lumerical simulations shown in Fig. 5 and Fig. 6 -the top and bottom layers were therefore completely misaligned. The results for these simulations are shown below in Fig. 7, where we define the polarization as TM for Fig. 7(a) and TE for Fig. 7(b).
Due to the average absorption of our device being above 80% (as given in Table 1 for all designs and polarizations), we choose this as the upper bound of the contour plot coloration for Figs. 7(a) and 7(b). Therefore, for absorption values 80% or above the color appears red, with diminishing absorption effects becoming more apparent with a change in color. Both Fig. 7(a) and 7(b) show very similar plots, with very broadband and highly-absorptive behavior occurring for both TM and TE polarizations over a large range of angles, where the incident angle was changed from 0°to 45°. Slightly better performance occurs for the TE polarized plot, which can be expected as the electric field of the incident light remains parallel to the resonant crosses and therefore allows correct operation with increasing angle of incidence -even at 45°-where the absorption bandwidth slightly broadens and blue-shifts, with a slight decrease in the absorption occurring at ∼2.75 THz. A similar case occurs for the TM polarized plot, with increasing oblique angle leading to further diminishing of the absorption around ∼2.75 THz, which expands to almost ∼4 THz at 45°. However, the TM polarized response is to red-shift with increasing angle rather than blue-shift (as for TE polarized). Interestingly, the bandwidth of the device seems to increase to the full width of the plot at 45°, corresponding to a bandwidth of over 6.5 THz, whilst still achieving greater than 50% absorption (yellow plot coloration). This lends itself to further showcase the benefits of our design, where neither the incident angle, polarization or layer alignment has too negative an effect on the absorption bandwidth or magnitude.

Conclusion
In summary we have detailed the design and simulation of significantly broadband THz absorbers formed of a bilayer fractal cross arrangement, which is not significantly affected by layer misalignments. We observe a record FWHM absorption of 104% (fractional bandwidth), whilst only being ∼λ/6 in thickness, where the average absorption (for four different designs) is ∼87%. Excellent performance is also shown for oblique angles up to 45°for both TE and TM polarizations, with a minimal decrease in the absorption efficiency and bandwidth. Taking into account the subwavelength thickness, layer-misalignment insensitivity, ultra-broadband and high absorption efficiency, and incident angle and polarization insensitivity, we propose that our work will be of benefit to realizing planar blackbody sources and absorbers, and of usefulness to developing absorber based sensing and imaging device. Due to the fabrication tolerances being more relaxed, as evidenced by the layer-misalignment insensitivity, the design shown will be of particular interest for commercialization.