Active spatial control of photothermal heating and thermo-actuated convective flow by engineering a plasmonic metasurface with heterodimer lattices

Zhimin Jing; Peihang Li; Cuiping Ma; Jiaying Wang; Roberto Caputo; Alexander O. Govorov; Arup Neogi; Hongxing Xu; Zhiming Wang

doi:10.1364/PRJ.472087

Abstract

Thermo-plasmonics, using plasmonic structures as heat sources, has been widely used in biomedical and microfluidic applications. However, a metasurface with single-element unit cells, considered as the sole heat source in a unit cell, functions at a fixed wavelength and has limited control over the thermo-plasmonically induced hydrodynamic effects. Plasmonic metasurfaces with metal disk heterodimer lattices can be viewed to possess two heat sources within a unit cell and are therefore designed to photo-actively control thermal distributions and fluid dynamics at the nanoscale. The locations of heat sources can be switched, and the direction of the convective flow in the central region of the unit cell can be reversed by shifting the wavelength of the excitation source without any change in the excitation direction or physical actuation of the structural elements. The temperature and velocity of a fluid are spatiotemporally controlled by the wavelength selectivity and polarization sensitivity of the plasmonic metasurface. Additionally, we investigate the effects of geometric parameters on the surface lattice resonances and their impact on the temperature and fluid velocity of the optofluidic system. Our results demonstrate excellent optical control of these plasmonic metasurface heating and thermal convection performances to design flexible platforms for microfluidics.

1. INTRODUCTION

Thermo-plasmonics opens up new opportunities for areas such as photothermal therapy [1,2], drug delivery [3], biosensing [4], solar-light harvesting [5], superhydrophobic coatings [6,7], trapping and rotating nanoparticles [8,9], high-resolution spectroscopies of nanoscale chirality [10,11], and manufacturing [12]. The photothermal effects in plasmonic structures generate micro- or nanoscale nonuniformity in thermal distribution, inducing hydrodynamic effects such as natural convection [13,14], Marangoni flow [15], thermophoresis [16], thermo-osmosis [17], thermoelectric pulling [18], and depletion forces [19]. However, controlling the photothermal energy distribution and fluid motion in systems based on plasmonic metasurface or photoactive material is challenging. Various means have been adopted to tackle this issue, such as modifying the spectral response from the surface plasmon resonance (SPR) band of plasmonic metasurfaces. For example, modifying the geometrical configuration (interparticle coupling) can change the near-field interaction of spatially localized electric fields, giving rise to an obvious shift of the SPR [20 –23]. However, this scheme usually needs controllable manipulation at nanoscale dimensions, which is passive and challenging to implement. In addition, the spectral response of plasmonic metasurfaces embedded within a photoactive or nonlinear medium such as liquid crystals, phase-change materials, and layered two-dimensional or III–V semiconductors can be modified by a change in their dielectric properties [24 –27]. Nevertheless, an external stimulus such as an electrical, magnetic, or electrochemical effect is necessary to induce a reversible change in the refractive index of the medium—a process that can be slow and inefficient compared to an optically induced effect.

Thermo-optical and light–matter interactions can be actively manipulated by modifying the localized SPR through an appropriate choice of the intrinsic property of the excitation source, such as polarization state, angle of incidence, phase, beam shape, or wavelength of the laser. The polarization state of incident light has been shown to modify the dielectrophoretic forces to manipulate biological cells [28]. The phase of incident coherent radiation was used to modulate electromagnetic force distribution and generate potential traps to manipulate nanoparticles [29]. A radially polarized beam can enhance the excitation of SPR and trapping efficiency [30]. However, the excitation of metal nanoparticles with characteristic localized SPR (LSPR) is restricted by the fixed resonant wavelength of the excitation source. The light absorption efficiency and, consequently, the localized heating are significantly lowered when the excitation source is detuned from this resonant LSP energy or directed away from the nanoparticles.

The range of SPR interactions has been enlarged by exploiting surface lattice resonances (SLRs)—a collective, coupled mode between the LSPR of the single nanostructures and the diffractive behavior of the periodic lattice [31 –34]. SLR structures are photonic lattices with narrow resonances that can be tuned to any desired operating wavelength and exhibit localized electric field and micrometer-scale delocalization [35,36]. The in-plane diffraction modes in a homogeneous environment or the evanescent diffraction orders termed Rayleigh anomalies (RAs) in an inhomogeneous environment result in an extended range of the excitation field of the SLR [37]. The orthogonal or parallel coupling within an SLR structure enables polarization control of incident light to manipulate the light–matter interaction within the meta-structure [38]. In this regard, our recent work proves that SLR has more advantages and application values than LSPR in optofluidics [39].

Sign up for Photonics Research TOC. Get the latest issue of Photonics Research delivered right to you！Sign up now

However, from the perspective of the morphology and lattice arrangement of metasurfaces, once the structures are designed, these SLR structures are passive because the lattice parameters of the plasmonic crystal will restrict the operational wavelength. And these structures with single-element unit cells are generally considered to be a heat source in a unit cell, operating with limited control of heat distributions and fluid dynamics. Also, the direction of the fluid flow is not easily manipulated flexibly in a symmetrical passive system. Recently, asymmetric nanostructures have been demonstrated to provide an additional degree of freedom to manipulate SLR and can be used as a sensing platform [40,41]. Asymmetric plasmonic heterostructures are also increasingly used to achieve nanoscale localized heat generation [42 –44]. Another recent approach uses superlattices, composed of finite arrays of nanoparticle patches grouped into microscale arrays to extend the scope from single-mode (SM) to multimode lattice resonances [45,46]. The presence of a multimode resonance with an asymmetric nanostructures array is likely to provide both wavelength and polarization sensitivity to SLR structures in addition to the usual wave vector sensitivity to incident radiation. Motivated by the possibility of realizing an optically active SLR structure for thermo-plasmonics, the design of a plasmonic metasurface composed of metal disk heterodimer (MDH) lattices is studied.

This work presents an active thermo-plasmonic device for wavelength- and polarization-controlled optical heating and thermal convection by exploiting the SLRs of MDH arrays. Specifically, spatially modulated nanoscale control of the photoinduced heating sites and thermally induced convection can be achieved by switching the incident wavelength to the respective SLR modes. The excitation wavelength can also modulate the convective flow rate of the fluid in the radial direction and its temperature in the orthogonal direction. The design of the SLR is investigated mathematically using the plasmon-hybridization technique. A parametric study of the various factors influencing the nanoscale temperature control and fluid flow in the MDH optofluidic system is presented. The proposed method for the active spatial control of photothermal heating and its induced hydrodynamic effects has a tremendous prospect in the application of optofluidic systems.

2. MODEL AND FORMALISM OF PHOTOTHERMAL HEATING AND THERMAL CONVECTION

The schematic of the proposed MDH arrays is shown in Fig. 1(a). A square array of MDH immersed in an index matching solution ( $n_{1} = 1.46$ ) is fabricated onto a ${SiO}_{2}$ substrate ( $n_{2} = 1.46$ ) to obtain a homogenous environment. The unit cell of the MDH arrays is shown in Figs. 1(b) and 1(c). The period of the square array is $P$ , the diameters of the two asymmetric disks (disk 1 and disk 2) are $d_{1}$ and $d_{2}$ , respectively, and the gap between disk 1 and disk 2 is $g$ . The thickness $h$ of each disk is 30 nm. As a constituent material, Au is chosen due to its chemical stability, relatively low loss in the near-infrared band, and ease of fabrication for realizing an actual optofluidic system [47,48]. The proposed MDH arrays can be fabricated using electron beam lithography, gold deposition, and liftoff processes.

Figure 1.(a) Schematic of metal disk heterodimer (MDH) arrays consisting of a square array of MDH deposited onto a glass substrate. The environment of the MDH is index-matched with the substrate using solution ( $n = 1.46$ ). (b), (c) Schematic of the unit cell of the periodic array with $d_{1} = 80 nm$ , $d_{2} = 120 nm$ , $g = 40 nm$ , $P = 450 nm$ , and $h = 30 nm$ , where $d$ is the diameter of the metal disk, 1 and 2 correspond to the two disks of MDH, $g$ is the gap between disk 1 and disk 2, $P$ is the period of the array, and $h$ is the height of the disk. The steam bubbles presented here are just to illustrate the fluid vividly; they have nothing to do with boiling the fluid.

Illuminated by a normally incident plane wave, extinction (sum of absorption and scattering) of MDH arrays is defined as 1-T, where T is the transmission through the MDH arrays. MDH arrays absorb incident light, resulting in heat dissipation, and the thermal generation is derived from the heat power density $q (r)$ inside MDH arrays: $q (r, t) = \frac{1}{2} Re [J^{*} (r, t) \cdot E (r, t)] = \frac{ω ε_{0}}{2} Im ε (ω) {| E_{ω} (r, t) |}^{2},$ (1)where $J^{*} (r, t)$ is the complex conjugate of the local induced current density, which generates heat losses via the Joule effect, $E (r, t) = Re [E_{ω} (r) \cdot e^{- i ω t}]$ is the electric field in the MDH arrays, $ε_{0}$ is vacuum permittivity, and $ε (ω)$ is the frequency-dependent relevant permittivity of the MDH arrays. Extensive studies proved that this heat generation in MDH arrays results in a temperature distribution $T (r, t)$ within the MDH arrays governed by the heat transfer equation [49,50] $ρ_{m} c_{m} \frac{\partial}{\partial t} T (r, t) - κ_{m} \nabla^{2} T (r, t) = q (r, t),$ (2)where $ρ_{m}$ , $c_{m}$ , $κ_{m}$ are density, specific heat capacity at constant pressure, and thermal conductivity of the metal, respectively, and $T (r, t)$ is the spatial temperature distribution. Baffou et al. pointed out that although $q (r)$ may be non-uniform depending on the shape of the nanostructure, the temperature distribution throughout the nanostructure is highly uniform due to the high thermal conductivity of the metal compared to the surrounding fluid [51]. Since light absorption outside the MDH arrays is almost negligible, the heat power density $q (r)$ equals zero, but convection may exist. At this point, the temperature distribution in the surrounding medium is governed by [49] $ρ_{s} c_{s} {\frac{\partial}{\partial t} T (r, t) + \nabla \cdot [T (r, t) v (r, t)]} - κ_{s} \nabla^{2} T (r, t) = 0,$ (3)where $ρ_{s}$ , $c_{s}$ , $κ_{s}$ are mass density, specific heat capacity at constant pressure, and thermal conductivity of the surrounding fluid, respectively, and $v (r, t)$ is the fluid velocity in the surrounding medium.

A change in temperature of the surrounding fluid will change its density. In general, the density will decrease as the temperature increases. As a result, the fluid at other locations will flow toward the center of mass reduction in a toroidal form due to the continuity of fluid, thus causing an upward generation of convection of the fluid. This fluid velocity profile is represented by the Navier–Stokes equation [14,52] $\frac{\partial}{\partial t} v (r, t) + [v (r, t) \cdot \nabla] v (r, t) = v_{s} \nabla^{2} v (r, t) + f_{t} [T (r, t)],$ (4)with continuity equation $\nabla \cdot v = 0$ , where $v_{s}$ and $f_{t}$ are, respectively, the kinematic viscosity of the fluid and the volume force (here only buoyancy) per unit mass due to the temperature gradient. The initial values of velocity are $v (r) = 0 nm / s$ . In nonisothermal flow calculation, the Boussinesq approximation is used to deal with the buoyancy force caused by the temperature difference to simplify the calculation [14]: $f_{t} (T) = g β_{s} [T (r, t) - T_{0}] \hat{z},$ (5)where $g$ , $β_{s}$ , $T_{0}$ , $\hat{z}$ are gravitational acceleration, dilatation coefficient of fluid, initial temperature, and upward $z$ -direction unit vector, respectively.

All calculations for the classical electromagnetic properties, temperature distributions, and fluid distributions have been carried out by commercial finite-element software (COMSOL Multiphysics v 5.5). By using modules of “electromagnetic waves, frequency domain,” “heat transfer in solids and fluids,” “laminar flow,” and Multiphysics coupled interfaces of “electromagnetic heating,” “nonisothermal coupling,” Eqs. (1 )–(5) are solved to investigate the multidisciplinary problem in the MDH arrays. $P$ , $d_{1}$ , $d_{2}$ , and $g$ are set to 450 nm, 80 nm, 120 nm, and 40 nm, respectively. MDH arrays are placed on a 1 μm thick substrate and immersed in 1 μm thick fluid. The incident plane wave is impinging from the fluid side to interact with the MDH arrays, exciting SLRs. The dielectric constants of Au are taken from the experimental data given by Johnson and Christy [53]. The details of the modeling and computational settings can be found in Appendix A.

3. RESULTS AND DISCUSSION

A. Excitation and Understanding of SLR

The excitation of the SLRs initiates the light-actuated modulation of the fluid for MDH arrays. Figure 2 depicts the optical properties of plasmonic lattices. The extinction spectrum of MDH arrays for the electric field parallel to the $y$ axis shows two obvious resonances at 679 nm (mode I) and 775 nm (mode II) [Fig. 2(a) red solid curve]. A different response is obtained when the electric field is parallel to the $x$ axis, showing just one obvious resonance at 786 nm [Fig. 2(a) black dashed line]. The calculated positions of these resonances both occur near the ( $\pm 1$ , 0) and (0, $\pm 1$ ) RAs at around 657 nm (blue dotted line), implying the presence of SLRs. Moreover, we can also observe the ( $\pm 1$ , $\pm 1$ ) RAs at around 464 nm (blue dotted line). Here, for a square lattice $P$ and homogeneous index $n$ , these positions of RAs can be obtained by [54] $λ_{i, j} = P n \frac{\sqrt{(i^{2} + j^{2}) - j^{2} \sin^{2} θ} \pm i \sin θ}{i^{2} + j^{2}},$ (6)where $θ$ is the angle of incidence, and $i$ and $j$ are integers corresponding to the diffracted order of the RAs. In addition, the extinction spectra of homodimer 1 (consisting of two identical disks 1) arrays and homodimer 2 (consisting of two identical disks 2) arrays are given in Fig. 2(b). One can see obvious resonance (at 685 nm) of homodimer 1 arrays (orange solid curve) and the other obvious resonance (at 768 nm) of homodimer 2 arrays (cyan solid curve) from Fig. 2(b). It is not difficult to find that the resonance wavelength separation of homodimers 1 and 2 arrays ( $Δ λ = 83 nm$ ) is smaller than MDH arrays ( $Δ λ^{'} = 96 nm$ ), which indicates that coupling behavior occurs in MDH. Figures 2(c)–2(f) show the electric field distributions of MDH at the two resonant wavelengths when the incident electric field is $y$ polarized. It can be seen that the electric field distributions vary with the incident wavelength. At 679 nm, the large electric field is distributed mainly around disk 1, while at 775 nm, the strong electric field is localized mostly around disk 2. The current density distributions (represented with white arrows) in the same direction on the two disks in the $x - y$ plane indicate that the two oscillations are in-phase at 679 nm [Fig. 2(c)]. This resonance is thereby considered to be the symmetric mode. At 775 nm, the current density distributions on the two disks are in opposite directions in the $x - y$ plane [Fig. 2(e)]. It implies that the two oscillations are out of phase, so the resonance can be considered the antisymmetric mode. In particular, some electric field distributions can be seen at the gap between the two disks, indicating the plasmon hybridization in MDH. Additionally, the localized field distributions on the heterodimer and delocalized field distributions along the $y$ direction in Figs. 2(c) and 2(e) reveal the coupling between LSP and ( $\pm 1$ , 0) RAs, thus confirming the generation of dual SLRs. Figures 2(g) and 2(h) show the electric field distributions of MDH for the $x$ -polarized electric field at resonance. It is observed that at 786 nm, most of the electric field is distributed on the surface of disk 2. Furthermore, current density distributions (represented with white arrows) on the two disks in the same direction indicate the bonding oscillations [Fig. 2(g)]. Similarly, the localized field distributions on the heterodimer and the delocalized field distributions along the $x$ direction are observed in Figs. 2(g) and 2(h). It indicates the coupling between LSP and ( $\pm 1$ , 0) RAs, thus confirming the generation of a single SLR.

Figure 2.Optical properties of plasmonic metasurface. (a) Extinction spectra of MDH arrays for $x$ - and $y$ -polarized incidence conditions. Inset: unit cell of MDH arrays. (b) Extinction spectra of homodimer 1 and homodimer 2 arrays for $y$ -polarized incidence conditions. Inset: unit cells of homodimer 1 and 2 arrays. (c)–(h) Electric field distributions of MDH at the two resonant wavelengths of $\sim 679 nm$ and $\sim 775 nm$ under $y$ polarization (c)–(f) and $\sim 786 nm$ under $x$ -polarization (g), (h). (i) Schematic plasmon hybridization of MDH.

Figure 2(i) depicts the schematic plasmon hybridization of MDH with different polarizations. For $y$ -polarized excitation, the vertical dipole of disk 1 and the quadrupole of disk 2 at 679 nm (investigated in Appendix B) indicate that the symmetric mode is derived from the plasmon hybridization between the dipolar resonance of disk 1 with the quadrupolar resonance of disk 2. Furthermore, its energy is closer to the dipolar resonance of disk 1. The vertical dipoles of disk 1 and disk 2 at 775 nm (investigated in Appendix B) indicate that the antisymmetric mode derives from plasmon hybridization between the dipolar resonance of disk 2 with the dipolar resonance of disk 1, and its energy is closer to the dipolar resonance of disk 2. For $x$ -polarized excitation, the horizontal dipoles of disk 1 and disk 2 at 786 nm (investigated in Appendix B) indicate that a bonding dipole–dipole mode appears. It results from the plasmon hybridization between the dipolar resonance of disk 1 with disk 2, and its energy is lower than the respective dipolar resonance of disk 1 and disk 2. The plasmon hybridization method is used to mathematically analyze the relationship between plasmon hybridization modes and the resonance of the two disks. For the $y$ -polarized excitation case, the symmetric mode $Φ_{SM}$ and antisymmetric mode $Φ_{AM}$ can be expressed by [55] $| Φ_{SM} ⟩ = τ | ϕ_{1} ⟩ + (\sqrt{1 - τ^{2}}) | ϕ_{2} ⟩,$ (7) $| Φ_{AM} ⟩ = - (\sqrt{1 - τ^{2}}) | ϕ_{1} ⟩ + τ | ϕ_{2} ⟩,$ (8)where $ϕ_{1}$ and $ϕ_{2}$ are plasmon modes of disk 1 and disk 2, and $τ$ is a constant that depends on the relative energies of plasmon modes for disk 1 and disk 2 and their interaction strength. Since $Φ_{SM}$ is closest in energy to $ϕ_{1}$ and $Φ_{AM}$ is closest in energy to $ϕ_{2}$ , it is not difficult to conclude that $ϕ_{1}$ contributes significantly more to $Φ_{SM}$ than $ϕ_{2}$ , and $ϕ_{2}$ contributes significantly more to $Φ_{AM}$ than $ϕ_{1}$ . For the $x$ -polarized excitation case, the bonding dipole–dipole mode $Φ_{BDD}$ is given by $| Φ_{BDD} ⟩ = α | ϕ_{1} ⟩ + (\sqrt{1 - α^{2}}) | ϕ_{2} ⟩,$ (9)where $α$ is a constant that depends on the relative energies of plasmon modes for disk 1 and disk 2 and their interaction strength. Note that $τ ≫ \sqrt{1 - τ^{2}}$ in Eqs. (7) and (8), and $α < \sqrt{1 - α^{2}}$ in Eq. (9).

B. Controlled Photothermal Heating and Convective Flow

Computational results in Fig. 3 provide an overview of the nanoscale spatial control of photothermal heating and thermo-actuated convection flow based on SLRs. Figure 3 explores the photoinduced control of temperature and thermally induced convection distributions in MDH by illuminating it with different polarizations and wavelengths of incident light. Based on typical values found in the literature [13,56], an incident light flux of $10^{4} W / {cm}^{2}$ is chosen. Halas and Nordlander et al. reported that the heat transfer process occurs in a very short period (on the order of ps-ns) [57], so to observe the heat transfer process more clearly, the data are processed for $t = 10 ns$ . The arrows on the streamlines denote the direction of motion of the fluid. Different colors represent different magnitudes of the surface temperature and fluid velocity. For $y$ -polarized incidence [Figs. 3(a)–3(d)], Figs. 3(a) and 3(b) show the horizontal ( $z = 30 nm$ ) and vertical ( $y = 0 nm$ and $x = 0 nm$ ) 2D slices of fluid convection overlaid on the temperature distributions at a wavelength of 679 nm. The magnitude of the temperature distributions shows a prominent feature; disk 1 has more pronounced photoinduced heat generation than disk 2. These inhomogeneous heat distributions at the nanoscale naturally induce asymmetric convective fluid distributions around the MDH. These effects are reversed in Figs. 3(c) and 3(d), where the horizontal ( $z = 30 nm$ ) and vertical ( $y = 0 nm$ and $x = 0 nm$ ) 2D slices of fluid convection overlaid on the temperature distributions at a wavelength of 775 nm are presented. In this case, a higher photoinduced temperature distribution of disk 2 than disk 1 can be observed. Irradiation at 775 nm also generates asymmetric convective fluid distributions around the MDH. Such obvious differences in temperature and flow field distributions at different incident wavelengths prove the wavelength-dependent photoinduced temperature and thermally induced convection at the nanoscale. These simulated observations indicate the possibility of switching the location of the heat source and thereby reversing the direction of the convective flow in a small area by tuning the excitation wavelength without any beam deflection or translation of the structural elements of the asymmetric structures. For the $x$ -polarized incidence case, Figs. 3(e) and 3(f) show the horizontal ( $z = 30 nm$ ) and vertical ( $y = 0 nm$ and $x = 0 nm$ ) 2D slices of fluid convection overlaid on the temperature distributions at a wavelength of 786 nm. One can observe nearly the same temperature and fluid flow distributions as in Figs. 3(c) and 3(d), but the amplitudes of temperature and velocity are smaller than in Figs. 3(c) and 3(d). The distinct differences in the temperature and flow field distributions at different incident polarizations also prove the polarization-dependent photoinduced temperature and thermally induced convection at the nanometer scale. Please note that the temperature and flow field distributions presented here are responses of a single cell, which we want to use to emphasize the response difference in more detail under different illumination cases. The differences in temperature and convection patterns between different wavelengths and polarizations (at the macroscopic scale) are also presented in Appendix E, where we approximately simulate the response of the MDH array system. These results also prove the wavelength- and polarization-dependent photoinduced temperature and thermally induced convection at the macroscopic scale within the excitation time of the order of ns.

Figure 3.Controlled photothermal heating and thermo-actuated convective flow of MDH under different incident wavelengths and polarizations at the nanoscale. (a)–(f) Photoinduced control of temperature and thermally induced convection distributions in MDH. (g) Velocity map of the fluid flow along the $x$ direction ( $y = 0 nm$ and $z = 30 nm$ ). Three different measuring points $P_{1} (0, 0, 30)$ , $P_{2} (- 200, 0, 30)$ , and $P_{3} (200, 0, 30)$ are marked in green solid circle, magenta solid triangle, and turquoise solid square, respectively. All data are given at the time $t = 10 ns$ for the incident light flux of $10^{4} W / {cm}^{2}$ .

To further elucidate photo-actively controlled fluid dynamics, Fig. 3(g) depicts the velocity map of the fluid along the $x$ direction ( $y = 0 nm$ and $z = 30 nm$ ) at different wavelengths and polarizations of incident light. The direction of fluid flow is indicated by positive and negative values, and the flow rate is represented by numerical values. The velocity components of three points $P_{1} (0, 0, 30)$ , $P_{2} (- 200, 0, 30)$ , and $P_{3} (200, 0, 30)$ along the $x$ direction are measured in different cases. One can observe that the velocity magnitude and direction of the three measuring points differ under different incident conditions. Interestingly, the direction of the velocity at $P_{1} (0, 0, 30)$ is reversed when the incident situation is changed; the magnitude of the velocity at $P_{2} (- 200, 0, 30)$ of 679 nm incidence under the $y$ polarization case is larger than the other two cases. However, the magnitude of the velocity at $P_{3} (200, 0, 30)$ of 775 nm incidence under the $y$ polarization case is greater than in the other two cases. Therefore, it can be concluded that when the excitation wavelength or polarization state changes, a reversal of the direction of fluid flow occurs at the location between the two heat sources, and switching occurs at the location of two heat sources. Therefore the wavelength and polarization state of incident light can be used as an additional degree of freedom to switch the location of the heat source and reverse the direction of the convective flow in part of the region.

To unveil the reasons behind this active spatial controllable phenomenon, we calculate the total power dissipation (i.e., total absorbed power) of the MDH in Fig. 4(a). Two distinct peaks occur at 679 nm and 775 nm due to the absorption of two disks at dual SLRs under $y$ -polarized illumination [Fig. 4(a) red solid curve]. Only one peak occurs at 786 nm due to the absorption of two disks at single SLRs under $x$ -polarized illumination [Fig. 4(a) black dashed line], indicating the total photothermal heating in the MDH. The contribution of the respective heat generation of the two disks to the total photothermal heating is estimated from the percentage of heat generation and heat distributions as calculated in Figs. 4(b)–4(d). In the case of illumination by $y$ -polarized light, Fig. 4(b) shows that from 678 nm to 710 nm, more than 90% of the total photothermal heating of the system is contributed by disk 1; however, when the wavelength is larger than 766 nm, more than 90% of the total heat of the system is contributed by disk 2. Correspondingly, as depicted in Fig. 4(d), at 679 nm, heat is predominantly concentrated on disk 1, while at 775 nm, heat is distributed mostly on disk 2. This wavelength-dependent heat generation explains the wavelength-dependent temperature and convection distributions in Fig. 3. In the case of $x$ -polarized illumination, Fig. 4(c) shows that when the wavelength is larger than 676 nm, more than 80% of the total heat of the system is contributed by disk 2. Correspondingly, as depicted in Fig. 4(d), at 786 nm, heat is concentrated mainly on the surface of disk 2. This polarization-dependent heat generation justifies the polarization-dependent temperature and convection distributions in Fig. 3.

Figure 4.Heat generation of MDH arrays. (a) Total absorbed power of MDH arrays under different polarizations as a function of wavelength. Percentage contribution of heat generation in disks 1 and 2 to the total heat generation of MDH arrays under $y$ -polarized illumination (b) and $x$ -polarized illumination (c). (d) Heat generation at wavelengths under different polarizations.

The double SLRs excited by $y$ -polarized incident light have more attractive application prospects than the single SLR excited by $x$ -polarized incident light. Therefore, our further studies are restricted to $y$ -polarized optical excitations.

C. Wavelength-Dependent Spatial–Temporal Temperature and Velocity Distributions

According to Eq. (1), the photoinduced heat generation depends on the resonance frequency. Therefore, we illuminate the MDH array with a $y$ -polarized plane wave at resonance frequencies around 679 nm and 775 nm to investigate the wavelength-dependent spatial–temporal temperature and velocity distributions of the MDH. First, Fig. 5 depicts the wavelength-dependent spatial properties. Figures 5(a) and 5(b) show schematics of $x - y$ and $x - z$ planes marked at three different locations. These planes pass through the center of disks 1 and 2, and between the two disks of the MDH, respectively. The parametric dependence of the spatial temperature distribution under different illumination wavelengths is shown in Figs. 5(c)–5(f). The results specifically show the incident wavelength-dependent spatial modulation of the nanoscale heating source. More specifically, at 679 nm, the temperature of disk 1 is much higher than that of disk 2. At 775 nm, the temperature of disk 1 is slightly lower than that of disk 2. These results confirm that dual SLRs can induce a spatial variation of temperature by active wavelength control and subsequently can induce controlled thermal convection in opposite directions.

Figure 5.Wavelength-dependent spatial axial temperature and velocity distributions in MDH. Schematic of the (a) $x - y$ plane and (b) $x - z$ plane marked with three lines of MDH. Temperature distributions of MDH under different incident wavelengths along the (c) $x$ direction and (d) $y$ direction. Temperature distributions of MDH along three different $z$ directions (lines 1, 2, and 3) at (e) 679 nm and (f) 775 nm. Velocity component distributions of MDH under different incident wavelengths along the (g) $x$ direction ( $y = 0 nm$ , $z = 458 nm$ ), (h) $y$ direction ( $x = 0 nm$ , $z = 458 nm$ ), and (i) $z$ direction (line 2). All data are acquired when the system reaches a steady state. The maximum temperatures along lines 1, 2, and 3 are 51.9°C, 45.5°C, and 41.1°C at 679 nm incident wave, and the maximum temperatures along lines 1, 2, and 3 are 39.4°C, 42.2°C, and 45.7°C at 775 nm incident wave, respectively.

Figure 5(c) shows that the maximum thermal energy is concentrated on disks 1 and 2 with a fairly uniform temperature distribution over the surface of the disks. This is due to the much higher thermal conductivity of Au with respect to its surrounding [58]. Similarly, from the wavelength-dependent spatial velocity distribution part [Figs. 5(g)–5(i)], one can see that the fluid convection appears at different velocities under different illumination wavelengths. Figure 5(i) plots the evolution of the fluid velocity component (measured at line 2) along the normal direction (as a function of the $z$ position). The velocity in the fluid along the $z$ direction exhibits an asymmetric Gaussian profile with a maximum velocity around $z = 458 nm$ . Please note that when the height of the fluid chamber changes, it still follows the same Gaussian trend, but the position of the maximum velocity will change accordingly with the change in the height of the fluid chamber. Figures 5(g) and 5(h) show the evolution of the fluid velocity component (measured at $z = 458 nm$ ) as a function of $x$ and $y$ positions, respectively. In Fig. 5(g), the fluid velocity changes from positive to negative values along the $x$ direction ( $y = 0$ , $z = 458 nm$ ), for cases of excitation light at both 679 nm and 775 nm. It is worth noting that the sections marked with dashed and solid purple lines indicate that the fluid velocity reverses between the case of excitation light at 679 nm and 775 nm. In Fig. 5(h), the fluid velocity changes from positive values to negative values along the $y$ direction ( $x = 0$ , $z = 458 nm$ ), for cases of excitation light at both 679 nm and 775 nm. The maximum velocity occurs in between the two disks. It is due to the heat transfer from the two asymmetric disks to the surrounding fluid. Since the maximum temperature (of disk 1) at the incident wavelength of 679 nm is higher than that (of disk 2) at 775 nm, the fluid velocity of the optofluidic system is also relatively faster at 679 nm irradiation than that at 775 nm.

The wavelength-dependent temporal temperature and velocity distributions of the MDH are shown in Fig. 6. The transient variation in the temperature of the two disks is shown in Figs. 6(a) and 6(b). The rate at which the fluid can be heated and optically driven is studied from the transient characteristics of the average temperature and the velocity of the fluid [Figs. 6(c) and 6(d)]. One can observe that the transient change in the velocity and temperature depends not only on the area of fluid over the nanoscale heater but also on the laser wavelength. Interestingly, with the increase of time, the temperature of the two disks and fluid gradually reaches a steady state at around 2 μs. The time-varying temperature behavior of fluid can be explained by [13] $T = T_{amb} + \frac{I A η}{B} + b e^{- B Δ t / c},$ (10)where $T_{amb}$ is environment temperature, $B$ is the thermal conductivity coefficient from fluid to the external environment, $I$ and $A$ are the incident light flux and the surface area exposed to radiation absorbance, respectively, $η$ is the photothermal conversion efficiency of the optofluidic system, $b$ is a fitting coefficient that can be obtained through specific experiment, and $Δ t$ and $c$ refer to the irradiation time and specific heat capacity of the solution, respectively. So when the irradiation time is long enough, there will be $e^{- B Δ t / c} \approx 0$ , and the temperature will reach a steady state. Additionally, the transient behavior of the velocity and temperature of the fluid evolves at a characteristic rate related to the steady-state profile of temperature and velocity in the same order as discussed in detail by Donner et al [14].

Figure 6.Wavelength-dependent temporal temperature and velocity distributions in MDH. (a) Temperature of two disks as a function of time at 679 nm. (b) Temperature of two disks as a function of time at 775 nm. (c) Average temperature and velocity of the fluid as a function of time at 679 nm. (d) Average temperature and velocity of the fluid as a function of time at 775 nm.

The spatial and temporal variations in temperature and velocity at the two SLR excitation wavelengths shown in Figs. 5 and 6 fully reflect the potential for an active modulation of the functionality of the MDH. Furthermore, the polarization state of incident light can be used to actively modulate the spatiotemporal change in temperature and velocity of the fluid, as shown in Appendix C. These figures show the fluid behavior for $x$ -polarized incident light. Thus, a novel optically driven non-contact mechanism for effectively heating and shifting the location of the heat source by modifying the color and polarization of the incident light source can be realized. It opens up a new method for developing flexible platforms in microfluidics applications.

D. Geometric Parametrization of the MDH-SLR

The lattice characteristics of the metasurface can modify the photothermal efficiency and the range of the thermal transport of fluid in a plasmonic SLR device. As with conventional SLRs, geometric parameters such as the lattice constant, scatterer dimension, heterodimer gap, and diameters of the heterodimer play a significant role in the thermal properties and fluid dynamics. Figures 7(a), 7(c), 7(e), and 7(g) show the extinction spectra as a function of the periodicity, gap, and diameters of the MDH. Figure 7(a) shows that as the period increases, the ( $\pm 1$ , 0), (0, $\pm 1$ ), and ( $\pm 1$ , $\pm 1$ ) RAs are excited and redshifted [the corresponding positions obtained from Eq. (6) are shown in dashed white and gray lines]. When the light is normally incident on the square lattice, ( $\pm 1$ , 0), (0, $\pm 1$ ) RAs are degenerate and coincide. The positions of the dual SLRs also exhibit a redshift due to the redshift of (0, $\pm 1$ ) and (0, $\pm 1$ ) RAs. Figure 7(c) shows that the dual SLR wavelengths are almost constant as the gap increases. It is attributed to an increase in the gap that does not substantially affect two oscillations along the $y$ direction. Moreover, a change in the isolated disk dimension will modify the symmetry of the MDH and thus change the extinction spectra. As revealed in Fig. 7(e), when the diameter of disk 1, $d_{1}$ , increases while disk 2, $d_{2}$ , remains unchanged at 120 nm, the SLR wavelength changes significantly. When $d_{1}$ is small, the intensity of mode I, which is dominated by disk 1, is also lower. As $d_{1}$ increases, mode I redshifts and increases in intensity. As it approaches the value of $d_{2}$ , the asymmetry is reduced, and mode II intensity gradually decreases. When $d_{1}$ equals $d_{2}$ , the heterodimer becomes a homodimer, resulting in the existence of an SM. Similarly, in Fig. 7(g), when the diameter of disk 2, $d_{2}$ , increases while disk 1, $d_{1}$ , remains unchanged at 80 nm, the SLR wavelength changes strikingly. When $d_{2}$ is similar to $d_{1}$ , the heterodimer becomes a homodimer, so only SM exists. With the gradual increase of $d_{2}$ , the asymmetry increases, and two SLR modes appear, which can be interpreted as mode splitting caused by plasmon hybridization. As $d_{2}$ continues to increase, mode II significantly redshifts with increased intensity. Besides geometric parameters, the effect of incident angles on the extinction spectra is also investigated in Appendix D.

Figure 7.Flexibility of resonance wavelength, temperature, and velocity. (a) Extinction spectra plotted as a function of wavelengths and periods. (b) Ratio of the temperature difference and fluid velocity of resonance mode I to mode II as a function of periods. (c) Extinction spectra plotted as a function of wavelengths and gaps. (d) Ratio of the temperature difference and fluid velocity of resonance mode I to mode II as a function of gaps. (e) Extinction spectra plotted as a function of wavelengths and diameters of disk 1. (f) Ratio of the temperature difference and fluid velocity of resonance mode I to mode II as a function of $d_{1}$ . (g) Extinction spectra plotted as a function of wavelengths and diameters of disk 2. (h) Ratio of the temperature difference and fluid velocity of resonance mode I to mode II as a function of $d_{2}$ .

The dual SLRs’ geometric parameters also regulate the temperature and fluid velocity. The ratio of the temperature difference and fluid velocity of resonance mode I to II has been plotted to demonstrate the controllability. Since the two disks have different surface temperatures under different resonance modes, we take only the absolute value of the temperature difference between the two disks in different resonance modes. As shown in Figs. 7(b), 7(d), 7(f), and 7(h), the ratios of the temperature difference and the fluid velocity are investigated. $Δ T$ represents the absolute value of temperature difference between the two disks, $V$ denotes the fluid velocity, and subscripts I and II represent resonance modes I and II. The value of $Δ T_{I} / Δ T_{II}$ indicates the ability of dual SLRs to adjust the surface temperature on the two disks (see orange lines and symbols). The value of $V_{I} / V_{II}$ shows the ability of dual SLRs to adjust the fluid velocity in the optofluidic system (see green lines and symbols).

Specifically, the value of $Δ T_{I} / Δ T_{II}$ decreases with the increase of period and gap [Figs. 7(b) and 7(d)]. It implies that mode I has a greater ability to regulate the surface temperatures of both disks than mode II at a smaller period and gap. At a larger period and gap, the situation is the opposite: mode II has a greater ability than mode I to regulate the surface temperatures of both disks. Moreover, it is observed from Fig. 7(b) that the ratio of temperature change ( $Δ T_{I} / Δ T_{II}$ ) follows the same trend as that of the change in velocity with the increasing lattice constant or periodicity of the plasmonic crystal. The difference in temperature and its effect on the fluid flow rate are not observed for periodicity less than 480 nm. However, the fluid flow rate is independent of the plasmon gap or the distance between the dimers, as shown in Fig. 7(d). It implies that the periodicity of the lattice has a stronger ability to regulate the fluid velocity than the gap. The dimension of the scatterer or the isolated disk diameter also affects $Δ T_{I} / Δ T_{II}$ and $V_{I} / V_{II}$ . As shown in Fig. 7(f), $Δ T_{I} / Δ T_{II}$ first increases and then rapidly decreases with the value of $d_{1}$ ; $V_{1} / V_{2}$ first increases slowly with $d_{1}$ and then increases rapidly when $d_{1}$ increases to 110 nm. These phenomena can be explained as follows: when $d_{1} = 60 nm$ , the intensity of mode I is much lower than that of mode II; thus, the values of $Δ T_{I} / Δ T_{II}$ and $V_{1} / V_{2}$ are both less than one, i.e., mode II has a stronger ability to regulate optical heating and thermal convection than mode I at $d_{1} = 60 nm$ . Since the absorption intensity increases with $d_{1}$ from 70 nm to 110 nm, mode I gradually becomes stronger in adjusting optical heating and thermal convection than mode II. A special case where double SLRs decompose into a single SLR occurs when $d_{1}$ equals $d_{2}$ , i.e., 120 nm. Figure 7(h) shows that $Δ T_{I} / Δ T_{II}$ and $V_{1} / V_{2}$ first decrease, then increase with the value of $d_{2}$ . It is noteworthy that the values of $Δ T_{I} / Δ T_{II}$ and $V_{1} / V_{2}$ are larger than one, indicating that mode I has a stronger ability to regulate optical heating and thermal convection than mode II when $d_{2}$ varies within a range from 80 nm to 160 nm. Also, a special case occurs when $d_{2}$ equals $d_{1}$ , i.e., 80 nm, where the double SLRs can be switched to a single SLR.

4. CONCLUSION

In conclusion, active and reversed spatial control of optical heating and thermal convection are realized on a plasmonic metasurface with heterodimer lattices by illuminating different light wavelengths or polarized states. Corresponding experimental possibilities are discussed in Appendix F. Thermodynamical and hydrodynamical results show that the different heat sources can be switched on and off. Thus inhomogeneous heat distributions at the nanoscale and spatial convective fluid distributions in a part of the region are obtained. Wavelength-dependent spatial and temporal temperature and velocity distributions of MDH are also investigated to demonstrate the controlled optical heating and thermal convection. Importantly, analytical and numerical results elucidate that the mechanism behind the regulation is the excitation of SLRs. The selectivity of the excitation wavelength can also control the spatiotemporal flow of the fluid. Additionally, the effects of geometric parameters on SLRs and corresponding photoinduced heat generation and thermally induced fluid motion are also analyzed. Our results highlight the controllable thermo-plasmonics potential of SLR-based metasurfaces and expand the capabilities of SLR-based metasurfaces for microfluidic applications.

APPENDIX A: SIMULATION SETTING METHOD

y - z

Figure 8. $x - z$ cross-sectional view of the model used in the simulation of the MDH.

Table 1.

Optical and Thermal Properties of Materials Used in the Calculations

Materials	Thermal Conductivity [W/(m·K)]	Mass Density [ $kg / m^{3}$ ]	Specific Heat Capacity [J/(kg·K)]	Refractive Index
Au	318	19320	129	Ref. [53]
Substrate (SiO2)	1.4	2650	840	1.46
Fluid (NaCl)	0.6	2003	1651	1.46

APPENDIX B: SURFACE CHARGE DISTRIBUTIONS OF MDH

To easily understand the underlying physics of the plasmon hybridization of MDH with different polarizations, surface charge distributions at the resonant wavelength under different polarization conditions are calculated, as shown in Fig. 9 . The red and blue parts in Fig. 9 denote the positive and negative charges, respectively. Surface charge distributions indicate that the vertical dipole of disk 1 and the quadrupole of disk 2 are excited at 679 nm [Fig. 9 (a)], the vertical dipoles of disk 1 and disk 2 are excited at 775 nm [Fig. 9 (b)], and the horizontal dipoles of disk 1 and disk 2 are excited at 786 nm [Fig. 9 (c)].

Figure 9.Surface charge distributions of resonance modes for different illuminated polarizations. Surface charge distributions of MDH at (a) 679 nm and (b) 775 nm for $y$ -polarized excitation case. (c) Surface charge distributions of MDH at 786 nm for $x$ -polarized excitation case.

APPENDIX C: SPATIOTEMPORAL TEMPERATURE AND VELOCITY DISTRIBUTIONS FOR X-POLARIZED INCIDENT LIGHT

x

Figure 10.Wavelength-dependent spatial axial temperature and velocity distributions in MDH in the case of $x$ -polarized incident wave. (a) Temperature distributions of MDH along the $x$ direction ( $y = 0 nm$ , $z = 0 nm$ ). (b) Velocity component distributions of MDH along the $x$ direction ( $y = 0 nm$ , $z = 458 nm$ ). (c) Temperature distributions of MDH along the $y$ direction ( $x = 0 nm$ , $z = 0 nm$ ). (d) Velocity component distributions of MDH along the $y$ direction ( $x = 0 nm$ , $z = 458 nm$ ). (e) Temperature distributions of MDH along different $z$ directions (lines 1, 2, and 3). (f) Velocity component distributions of MDH along the $z$ direction (line 2). Schematics of the $x - y$ plane and $x - z$ plane marked with three lines of MDH are shown in Figs. 5(a) and 5(b). The maximum temperatures along lines 1, 2, and 3 are 38.0°C, 39.6°C, and 42.0°C at 786 nm incident wave, respectively.

Figure 11.Wavelength-dependent temporal temperature and velocity distributions in MDH in the case of $x$ -polarized incident wave. (a) Temperature of two disks as a function of time under the incident wavelength of 786 nm. (b) Average temperature and velocity of fluid as a function of time under the incident wavelength of 786 nm.

APPENDIX D: SPATIOTEMPORAL TEMPERATURE AND VELOCITY DISTRIBUTIONS

Besides geometric parameters, the effect of incident angles on the extinction spectra is also investigated. Figure 12 shows the angle at which the transition from SM to multimode SLR can be realized. Thus, the momentum conservation of the MDH-SLR structure can be applied for regulating temperature and fluid velocity and opening up routes for the multimode operation of nanolasers.

Figure 12.Extinction spectra of the MDH with different incident angles of $α = 0 °$ , 2°, 4°, 6°, 8°, 10°. With the increase of $α$ , mode II exhibits an obvious redshift, mode I gradually splits, and more new modes appear. It is conceivable that this phenomenon will inevitably lead to the flexibility of optical heating and thermal convection.

APPENDIX E: PHOTOTHERMAL HEATING AND THERMAL CONVECTION RESPONSE AT MACROSCOPIC SCALE

5 \times 5

Figure 13.Temperature and fluid convection patterns under different incident wavelengths and polarizations of $5 \times 5$ array of MDH. The data are given at different illumination times ( $t = 10 ns$ , 1000 ns).

APPENDIX F: EXPERIMENTAL POSSIBILITIES

The experimental setup for thermal convection flow refers to that used in some previous studies [60–62]. The thermo-plasmonic device can be driven using a diffraction-limited or super-resolution microscope. The motion of the fluid can be simultaneously observed under white light illumination. The collimated or focused beam from the microscope incident on the heterodimer disk can be used to generate the convection flow. A fast dynamic CMOS camera can be used to record the fluid flow and the convection pattern.

x - y

In terms of illumination time, as the time changes from the order of ns to the order of μs, the system gradually tends to a steady state. The system will eventually display similar Raleigh-Benard-like fluid convection, but with different magnitudes of temperature, under different excitation cases. Given this, we can use a periodic nanosecond pulse laser as the fluid pump source to achieve wavelength- and polarization-controlled temperature and flow field distribution. In general, the demonstration and quantification of this thermally actuated convection provide the feasibility for optical control of the plasmon-assisted experiment, such as driving and trapping.

References

[1] Y. Liu, F. Mo, J. Hu, Q. Jiang, X. Wang, Z. Zou, X.-Z. Zhang, D.-W. Pang, X. Liu. Precision photothermal therapy and photoacoustic imaging by in situ activatable thermoplasmonics. Chem. Sci., 12, 10097-10105(2021).

[2] X. Huang, F. Zha, J. Zou, Y. Li, F. Wang, X. Chen. Photoacoustic imaging-guided synergistic photothermal/radiotherapy using plasmonic Bi/Bi₂O_3−x nanoparticles. Adv. Funct. Mater., 32, 2113353(2022).

[3] R. Huschka, A. Barhoumi, Q. Liu, J. A. Roth, L. Ji, N. J. Halas. Gene silencing by gold nanoshell-mediated delivery and laser-triggered release of antisense oligonucleotide and siRNA. ACS Nano, 6, 7681-7691(2012).

[4] G. Qiu, Z. Gai, Y. Tao, J. Schmitt, G. A. Kullak-Ublick, J. Wang. Dual-functional plasmonic photothermal biosensors for highly accurate severe acute respiratory syndrome coronavirus 2 detection. ACS Nano, 14, 5268-5277(2020).

[5] M. Zhu, Y. Li, F. Chen, X. Zhu, J. Dai, Y. Li, Z. Yang, X. Yan, J. Song, Y. Wang, E. Hitz, W. Luo, M. Lu, B. Yang, L. Hu. Plasmonic wood for high-efficiency solar steam generation. Adv. Energy Mater., 8, 1701028(2018).

[6] B. Wang, Z. Jing, M. Zhao, P. Yu, E. Ashalley, P. Li, C. Ma, X. Tong, R. Caputo, A. O. Govorov, H. Xu, Z. M. Wang. Ultraflexible photothermal superhydrophobic coating with multifunctional applications based on plasmonic TiN nanoparticles. Adv. Opt. Mater., 10, 2200168(2022).

[7] B. Wang, P. Yu, Q. Yang, Z. Jing, W. Wang, P. Li, X. Tong, F. Lin, D. Wang, G. E. Lio, R. Caputo, O. Ávalos-Ovando, A. O. Govorov, H. Xu, Z. M. Wang. Upcycling of biomass waste into photothermal superhydrophobic coating for efficient anti-icing and deicing. Mater. Today Phys., 24, 100683(2022).

[8] X. Zhao, Y. Shi, T. Pan, D. Lu, J. Xiong, B. Li, H. Xin. In situ single-cell surgery and intracellular organelle manipulation via thermoplasmonics combined optical trapping. Nano Lett., 22, 402-410(2022).

[9] L. Shao, M. Käll. Light-driven rotation of plasmonic nanomotors. Adv. Funct. Mater., 28, 1706272(2018).

[10] Y. Chen, W. Du, Q. Zhang, O. Ávalos-Ovando, J. Wu, Q.-H. Xu, N. Liu, H. Okamoto, A. O. Govorov, Q. Xiong, C.-W. Qiu. Multidimensional nanoscopic chiroptics. Nat. Rev. Phys., 4, 113-124(2022).

[11] O. Ávalos-Ovando, E. Y. Santiago, A. Movsesyan, X.-T. Kong, P. Yu, L. V. Besteiro, L. K. Khorashad, H. Okamoto, J. M. Slocik, M. A. Correa-Duarte, M. Comesaña-Hermo, T. Liedl, Z. Wang, G. Markovich, S. Burger, A. O. Govorov. Chiral bioinspired plasmonics: a paradigm shift for optical activity and photochemistry. ACS Photon., 9, 2219-2236(2022).

[12] A. W. Powell, A. Stavrinadis, I. de Miguel, G. Konstantatos, R. Quidant. White and brightly colored 3D printing based on resonant photothermal sensitizers. Nano Lett., 18, 6660-6664(2018).

[13] C. Ma, P. Yu, W. Wang, Y. Zhu, F. Lin, J. Wang, Z. Jing, X.-T. Kong, P. Li, A. O. Govorov, D. Liu, H. Xu, Z. Wang. Chiral optofluidics with a plasmonic metasurface using the photothermal effect. ACS Nano, 15, 16357-16367(2021).

[14] J. S. Donner, G. Baffou, D. McCloskey, R. Quidant. Plasmon-assisted optofluidics. ACS Nano, 5, 5457-5462(2011).

[15] F. Winterer, C. M. Maier, C. Pernpeintner, T. Lohmüller. Optofluidic transport and manipulation of plasmonic nanoparticles by thermocapillary convection. Soft Matter, 14, 628-634(2018).

[16] S. Duhr, D. Braun. Optothermal molecule trapping by opposing fluid flow with thermophoretic drift. Phys. Rev. Lett., 97, 038103(2006).

[17] A. P. Bregulla, A. Würger, K. Günther, M. Mertig, F. Cichos. Thermo-osmotic flow in thin films. Phys. Rev. Lett., 116, 188303(2016).

[18] L. Lin, P. S. Kollipara, A. Kotnala, T. Jiang, Y. Liu, X. Peng, B. A. Korgel, Y. Zheng. Opto-thermoelectric pulling of light-absorbing particles. Light Sci. Appl., 9, 34(2020).

[19] P. Liu, S. Ye, F. Ye, K. Chen, M. Yang. Constraint dependence of active depletion forces on passive particles. Phys. Rev. Lett., 124, 158001(2020).

[20] O. A. Yeshchenko, V. V. Kozachenko, A. V. Tomchuk, M. Haftel, R. J. Knize, A. O. Pinchuk. Plasmonic metasurfaces with tunable gap and collective surface plasmon resonance modes. J. Phys. Chem. C, 123, 13057-13062(2019).

[21] K. Martens, F. Binkowski, L. Nguyen, L. Hu, A. O. Govorov, S. Burger, T. Liedl. Long- and short-ranged chiral interactions in DNA-assembled plasmonic chains. Nat. Commun., 12, 2025(2021).

[22] A. Ferraro, G. E. Lio, A. Hmina, G. Palermo, J. M. Djouda, T. Maurer, R. Caputo. Tailoring of plasmonic functionalized metastructures to enhance local heating release. Nanophotonics, 10, 3907-3916(2021).

[23] E. Ashalley, C.-P. Ma, Y.-S. Zhu, H.-X. Xu, P. Yu, Z.-M. Wang. Recent progress in chiral absorptive metamaterials. J. Electron Sci. Technol., 19, 100098(2021).

[24] N. H. Tu, K. Yoshioka, S. Sasaki, M. Takamura, K. Muraki, N. Kumada. Active spatial control of terahertz plasmons in graphene. Commun. Mater., 1, 7(2020).

[25] Y. Wang, D. Wei, P. Sohr, J. M. O. Zide, S. Law. Extending the tunable plasma wavelength in III–V semiconductors from the mid-infrared to the short-wave infrared by embedding self-assembled ErAs nanostructures in GaAs. Adv. Opt. Mater., 8, 1900937(2020).

[26] V. K. S. Hsiao, Y. B. Zheng, B. K. Juluri, T. J. Huang. Light-driven plasmonic switches based on Au nanodisk arrays and photoresponsive liquid crystals. Adv. Mater., 20, 3528-3532(2008).

[27] S. Rubin, Y. Fainman. Nonlinear, tunable, and active optical metasurface with liquid film. Adv. Photon., 1, 066003(2019).

[28] X. Miao, L. Y. Lin. New opto-plasmonic tweezers for manipulation and rotation of biological cells-design and fabrication. International Conference of the IEEE Engineering in Medicine and Biology Society, 4318-4321(2006).

[29] Y. Zhang, W. Shi, Z. Shen, Z. Man, C. Min, J. Shen, S. Zhu, H. P. Urbach, X. Yuan. A plasmonic spanner for metal particle manipulation. Sci. Rep., 5, 15446(2015).

[30] X. Wang, Y. Zhang, Y. Dai, C. Min, X. Yuan. Enhancing plasmonic trapping with a perfect radially polarized beam. Photon. Res., 6, 847-852(2018).

[31] Z. Li, S. Butun, K. Aydin. Ultranarrow band absorbers based on surface lattice resonances in nanostructured metal surfaces. ACS Nano, 8, 8242-8248(2014).

[32] M. S. Bin-Alam, O. Reshef, Y. Mamchur, M. Z. Alam, G. Carlow, J. Upham, B. T. Sullivan, J.-M. Ménard, M. J. Huttunen, R. W. Boyd, K. Dolgaleva. Ultra-high-Q resonances in plasmonic metasurfaces. Nat. Commun., 12, 974(2021).

[33] A. K. Boddeti, J. Guan, T. Sentz, X. Juarez, W. Newman, C. Cortes, T. W. Odom, Z. Jacob. Long-range dipole–dipole interactions in a plasmonic lattice. Nano Lett., 22, 22-28(2022).

[34] Y. Lin, D. Wang, J. Hu, J. Liu, W. Wang, J. Guan, R. D. Schaller, T. W. Odom. Engineering symmetry-breaking nanocrescent arrays for nanolasing. Adv. Funct. Mater., 29, 1904157(2019).

[35] B. B. Rajeeva, L. Lin, Y. Zheng. Design and applications of lattice plasmon resonances. Nano Res., 11, 4423-4440(2018).

[36] V. G. Kravets, A. V. Kabashin, W. L. Barnes, A. N. Grigorenko. Plasmonic surface lattice resonances: a review of properties and applications. Chem. Rev., 118, 5912-5951(2018).

[37] X. Yang, G. Xiao, Y. Lu, G. Li. Narrow plasmonic surface lattice resonances with preference to asymmetric dielectric environment. Opt. Express, 27, 25384-25394(2019).

[38] L. Lin, Y. Zheng. Engineering of parallel plasmonic–photonic interactions for on-chip refractive index sensors. Nanoscale, 7, 12205-12214(2015).

[39] Z. Jing, P. Yu, A. Movsesyan, C. Ma, P. Li, Y. Zhu, A. O. Govorov, A. Neogi, Z. Wang. Manipulation of fluid convection by surface lattice resonance. Adv. Opt. Mater., 10, 2201066(2022).

[40] B. D. Thackray, V. G. Kravets, F. Schedin, G. Auton, P. A. Thomas, A. N. Grigorenko. Narrow collective plasmon resonances in nanostructure arrays observed at normal light incidence for simplified sensing in asymmetric air and water environments. ACS Photon., 1, 1116-1126(2014).