Effective surface conductivity of optical hyperbolic metasurfaces: from far-field characterization to surface wave analysis

Metasurfaces offer great potential to control near- and far-fields through engineering optical properties of elementary cells or meta-atoms. Such perspective opens a route to efficient manipulation of the optical signals both at nanoscale and in photonics applications. In this paper we show that a local surface conductivity tensor well describes optical properties of a resonant plasmonic hyperbolic metasurface both in the far-field and in the near-field regimes, where spatial dispersion usually plays a crucial role. We retrieve the effective surface conductivity tensor from the comparative analysis of experimental and numerical reflectance spectra of a metasurface composed of elliptical gold nanoparticles. Afterwards, the restored conductivities are validated by semi-analytic parameters obtained with the nonlocal discrete dipole model with and without interaction contribution between meta-atoms. The effective parameters are further used for the dispersion analysis of surface plasmons localized at the metasurface. The obtained effective conductivity describes correctly the dispersion law of both quasi-TE and quasi-TM plasmons in a wide range of optical frequencies as well as the peculiarities of their propagation regimes, in particular, topological transition from the elliptical to hyperbolic regime with eligible accuracy. The analysis in question offers a simple practical way to describe properties of metasurfaces including ones in the near-field zone with effective conductivity tensor extracting from the convenient far-field characterization.


S1. Distribution of nanodisks sizes
The fabrication of a metasurface is still challenged and complicated technological process. Obviously, not all particles have identical parameters, i.e. there is a distribution of particles position and sizes. According to such a distribution shown in

S2. Polarizability of thin nanodisk with elliptical base
We define the polarizability of a nanodisk with elliptical baseα 0 = diag(α x , α y , α z ) through the polarizability of an ellipsoid. First, we consider the case of an ellipsoid with semiaxes b x , b y , b z . For an anisotropic particle the depolarization factor N i should be introduced 1 : , i = x, y, z. (S1)

1/7
Three depolarization factors for any ellipsoid satisfy the relation: N x + N y + N z = 1. Finally, the polarizability of the ellipsoid is Here ε m is the permittivity of the surrounding medium and ε(ω) is the permittivity of the scatterer material. In case of a sphere, when b x = b y = b z , the depolarization factor is N i = 1/3 and we get the polarizability of a sphere according to Clausius-Mossotti relation 2 . Then, we switch from the ellipsoid to the elliptical nanodisk with sizes a x , a y , a z = H/2 and make the substitution for the semiaxes b i = (1.5) 1/3 a i , that takes into account the difference between the volumes of an ellipsoid and an elliptical cylinder. After that we use Eq. (S2) as the polarizability of the elliptical cylinder. Thus, we obtain the polarizability of a thin nanodisk with elliptical base (Fig. S2). We conclude that resonance of the normal component of polarizability α z is on very small wavelengths and for the studied range α z is much smaller that inplane components of polarizability. So, we can consider effective polarizability of thin nanodisk with elliptical base as a two-dimensional polarizability or conductivity tensor (5).

S3. Extraction of surface conductivity dispersion using fitting with Drude-Lorentz formula
We introduce the Drude-Lorentz model with three undefined coefficients β , γ, σ 0 : where σ 0 is the amplitude of conductivity dispersion, β = ω − Ω 2 /ω, ω is the operating frequency, Ω is the spectral position of resonance, γ is the bandwidth of resonance. Then we use Fresnel equations for a two-dimensional layer with an effective conductivity σ located between two isotropic media n 1 and n 2 in order to express the reflection coefficient: Here we use Gauss units and express surface conductivity in the normalized dimensionless units σ = 4π σ /c. Substituting Eq. (S3) into Eq. (S4) we obtain or expressing it through the reflectance We consider a metasurface under consideration ( Fig. 1) with the corresponding reflectance dispersion (Fig. 2a). Then we perform the fitting of the Eq. (S6) based on the least-squares method with Eq. (

3/7
the reflectance intensity in the simulation we increase the losses of the gold by six times (Fig. S3b). Hereinafter, we use the model of the gold permittivity from Ref. [3]. Finally, we obtain the dispersion of the surface conductivity for both polarizations from the fitting of the experimental data according to the Eq. (S3) [See Fig. S3c]. By knowing the parameters of fitting formula Eq. (S3) we can explicitly find the phase of the reflection coefficient from Eq. (S4) [See Fig. S3a].

S4. Performing surface conductivity extraction with oblique angle of incidence
The retrieval procedure of the effective surface conductivity within zero-thickness approximation can be applied for the arbitrary angle of incidence by modifying the Fresnel equations appropriately. Importantly, the final equations for the surface conductivity split by polarizations now and depend on the incident θ i and refraction θ t angles 4 : σ y = n 1 cos θ i − n 2 cos θ t − S y 11 (n 1 cos θ i − n 2 cos θ t ) 1 + S y 11 (S7) These formulas are highly important in the case of the strong spatial dispersion, because they allow to determine the dependence of the retrieved conductivity on the in-plane wavevector.

S5. Phase correction in frame of zero-thickness approximation
Within zero-thickness approximation we substitute a plasmonic resonant metasurface of finite thickness H by a two-dimensional layer with effective conductivity σ . The obvious question arises: at which distance H eff from the substrate we should dispose a two-dimensional layer? To define this distance we use a single criterion connected to the energy conservation law: The total phase of reflection coefficient (S-parameter) obtained by simulation in CST Microwave Studio is composed of two terms. The first one is an intrinsic phase associated directly with the reflection from a metasurface ϕ 0 , while the second term is an extrinsic part caused by the electromagnetic waves propagation from the excitation port to metasurface and back (L 0 is a distance between port and top of a metasurface). It is extremely important to define the extrinsic phase and make the appropriate phase correction in order to obtain the reflection coefficient related to the metasurface properties intrinsically. So, we can express the total phase ϕ as where k 0 = nω/c, n is a refractive index of the super-or substrate. Figure S4 shows the real parts of conductivity tensor components for the different locations of a two-dimensional layer (top, middle and bottom of a metasurface of finite thickness). One can see that this criterion is satisfied only, when a two-dimensional layer is disposed at the distance H/2 from the substrate.

S6. Ewald summation for 2D square lattice
The lattice sums can be evaluated by using scalar Green's function in r-space and k-space: (S10) Here k 0 is the wavevector of light in free space, b is the reciprocal wavevector of the structure, r is the excitation dipole or structural defect position and r i is the position of i-th dipole of lattice structure. We consider the case r = (0, 0, 0). Scalar Green's function (S10) can be divided into two parts by using Ewald summation with Ewald parameter K: For two-dimensional layer in the xy-plane with square lattice a the first term is defined through the function where while the second term is expressed with the function The overall sum should be not significantly dependent on the Ewald parameter K. This parameter is taken as K ∼ 1/a, where a is the lattice constant. So, we represent the Green's function as a sum of two contributions. The first term is calculated in real space, while the second is calculated in k-space using Fourier transform. It significantly reduces calculation time, keeping accuracy up to 10 −4 . 5

S7. Impact of the difference between substrate and superstrate permittivities
The numerical and theoretical calculations have been performed under assumption of the homogeneous ambient medium. Indeed, the permittivities of the fused silica substrate and the immersion liquid may be slightly different in real experiment. Here, we analyze the impact of the perturbatively small difference between the permittivities of the substrate ε 0 and superstratẽ ε = ε 0 (1 + δ ε). Namely, we consider several cases of the permittivity difference from δ ε = −0.2 to δ ε = 0.1. The permittivity of substrate is chosen as the permittivity of the fused silica ε 0 = 2.1.
The simulation in CST Microwave Studio of the reflectance under normal incidence for both polarizations with different δ ε is shown in Figs. S5a,b. One can see that the bigger superstrate permittivity over substrate one leads to the redshift of the spectrum, while the blueshift takes place for the opposite situation. The dependences of the reflectance resonance position on the permittivity perturbation are presented in Fig. S5c for TM polarization and Fig. S5d for TE polarization. The dependence of the resonance position on the difference of the permittivities is close to linear. The exact coincidence of the resonances obtained experimentally and numerically has been achieved for δ ε = −0.2 in TM polarization and δ ε = −0.16 in TE polarization.
The modification of the discrete dipole model on the different substrate and superstrate for the surface conductivity tensor can be fulfilled by using the classical perturbation theory assuming δ ε 1. It results in the additional term in the conductivity tensor: The correction factor of the permittivity perturbation for the effective surface polarizability δα eff can be expressed through the corresponding perturbative terms for the disk polarizability and the interaction constant (as the substrate contribution to the Green's function) as follows: (S16) Here,Î is the two-dimensional unit dyad, δ jk is the Kronecker symbol. 6/7 Figure S6 shows the effective surface conductivity spectra for the fairly different substrate and superstrate materials (δ ε = 0.05) calculated numerically and by using the DDM with perturbation corrections. One can see that the obtained results are in good agreement.
To conclude, we have shown that the effect of the difference between substrate and immersion liquid permittivities does not have a significant impact on the discussed results, but can lead to the better correspondence between the numerical calculations and the experimental measurements.