Optical forces and trapping potentials of a dual-waveguide trap based on multimode solid-core waveguides

M. M. van Leest m.m.vanleest@tudelft.nl Delft University of Technology, Kavli Institute of Nanoscience, Department of Quantum Nanoscience, Lorentzweg 1, 2628 CJ Delft, The Netherlands F. Bernal Arango Delft University of Technology, Kavli Institute of Nanoscience, Department of Quantum Nanoscience, Lorentzweg 1, 2628 CJ Delft, The Netherlands Currently with FOM Institute AMOLF, P.O. Box 41883, 1009 Amsterdam, The Netherlands J. Caro Delft University of Technology, Kavli Institute of Nanoscience, Department of Quantum Nanoscience, Lorentzweg 1, 2628 CJ Delft, The Netherlands


INTRODUCTION
The combination of photonics with a lab-on-a-chip attracts strong attention in view of its applicability in bio-sensing [1].For example, combining laser tweezing and Raman spectroscopy [2], enabling both trapping and identification of single cells in a fluid, is an important development.New in this is miniaturization of the tweezers and their integration in the fluidic chip using optical waveguides.This leads to an integrated photonics version of Ashkin's dual-beam trap [3] which avoids cumbersome free-space laser beams and lenses.In dual-beam traps the scattering force acts in the axial direction of the beams and the gradient force in the transverse direction [3].Various dual-waveguide traps have been integrated in a fluidic chip, either using fibers [4]- [6] or liquid core waveguides [7].Although these prototype dual-waveguide traps are important proofs-of-principle of the integrated concept, they do not permit high volume production with standard lab-on-a-chip technology, a requirement for widespread application in areas such as water-quality monitoring and medical diagnostics.
Here, we propose a novel design of the dual-waveguide trap.Contrary to other approaches, we base the trap on microfabricated embedded solid-core waveguides.The waveguide material is Si 3 N 4 , while the cladding is SiO 2 .These choices enable fabrication fully compatible with mass-production of standard glass-based microfluidics, ensuring reproducible quality and robustness of the devices and enabling widespread application.
To study our design, we perform three-dimensional (3D) finite-difference time-domain (FDTD) simulations of the beam exiting a waveguide.Further, for the complete trap and for polystyrene beads, we calculate optical forces and trapping potentials.This is done for the Raman excitation wavelength λ = 785 nm, so that trap characteristics are evaluated under conditions for which also a single particle Raman spectrum can be measured with a standard spectrometer (see Section 6).

APPROACH
FIG. 1 Sketch of the dual-beam trap based on Si 3 N 4 solid-core waveguides.The origin of the coordinate system is in the trap center.
Figure 1 depicts the trap geometry, showing waveguides of thickness t and width w, and a gap L between the facets.Si 3 N 4 is chosen as waveguide material as it is very suitable for near-visible and visible Raman excitation wavelengths and for its compatibility with glass-based microfluidic technology.Furthermore, Si 3 N 4 has a good index contrast with SiO 2 (n Si3N4 = 1.9 versus n SiO2 = 1.45), which is the cladding material (see Figure 1).The SiO 2 layer has a trench of which the walls are in the same vertical plane as the waveguides facets, thus defining the trapping geometry in the fluidic channel.The width of the channel (equal to gap L) is chosen to be 5 µm.This size is compatible with various types op biological objects, for instance certain bacteria.It is noted that proper trap operation is not limited to this gap size.

PROFILE OF THE OUTPUT BEAM
The remaining variables are the waveguide thickness and width.Many biological objects such as bacteria are close to spherical, suggesting waveguides giving a highly symmetric trapping potential.Thus, we adopt a square cross section for the waveguide.We choose the side of the square to be 1 µm, so that the waveguide can accommodate up to 8 modes per polarization direction for λ = 785 nm [8].The choice of the waveguide cross section is discussed below, in connection to Figure 2.
As the mode structure of the waveguide cannot be derived analytically [9], we perform 3D FDTD simulations with Meep [10].The continuous excitation source is y-polarized.Its position is 4 µm from the facet and it delivers a power of P = 0.5 W. The medium in the fluidic channel is water, with n H2O = 1.33.The size of the calculation cell in the x-, yand zdirection is 17 µm, 11 µm and 11 µm, respectively.We use 2 µm thick perfectly matched layers.
Figure 2 shows characterics of the simulated beam exiting a single waveguide (taking the left beam in Figure 1, thus fixing it with respect to the axes).In Figure 2(a) the red color represents a time-averaged intensity I of the electrical field E(I ∝ E 2 )in the xy-plane.Most important in Figure 2(a) is the beam spreading, resulting from diffraction of modes at the waveguide facet and leading to a non-zero scattering force of the trap [3].Further, owing to beam spreading the gradient force acts in a wide region, enabling trapping of larger objects.The pattern in the waveguide shows a short period modulation on top of a wide maximum (see also Figure 2(c)).The modulation results from interference of the wave coming from the source and the wave partially reflected at the facet.The wide maximum is part of the multimode interference (MMI) pattern [11] occurring as a result of excitation of several modes in the waveguide.From other simulations we find that smaller cross sections give stronger beam spreading, while larger cross sections give MMI patterns outside the waveguide with a multiple ray structure, each situation preventing proper trapping characteristics.Among the sizes simulated, a square cross section of 1×1 µm 2 gives the optimum result.
In Figures 2(b) and 2(c) the averaged profiles on the y-axis and the x-axis, respectively, show a smooth bell-shaped curve and a smoothly decaying curve.However, before it decays, in the x-direction the beam shows a maximum, which is the continuation of the MMI pattern in the waveguide.Details of these curves depend on the source-facet distance, which influences the MMI pattern.The main characteristics leading to a well defined optical trap, viz. the bell shape and the decay, are seen for each source position.
Next, the averaged intensity distribution is simulated for the complete trap, with a 1 µm diameter polystyrene bead centered at the point (0, 1, 0), again using a 4 µm sourcefacet distance.The sources each deliver P = 0.5 W to a waveguide (i.e. total input power is 1 W) and they oscillate in phase.The results are shown in Figure 3, for the yz-plane (a) and for the xy-plane (b).The period of the intensity modulation in the water is 297 nm, in agreement with λ/2 = 785 nm/(2×1.33)= 295 nm.The effect of the bead on the distribution is apparent.As a result of the higher index of the bead (n polyst.= 1.6), it concentrates the field inside its volume, thus perturbing the distribution in its vicinity.This is also re-flected in the shorter period of the pattern in the bead.The global intensity distribution does not depend on the relative phase, as we find from simulations for various phase differences between the sources.The positioning of the periodicity in the gap, however, does depend on the relative phase.We also perform simulations for a z-polarized source.The resulting intensity distributions are very close to those of Figures 2  and 3, as expected on the basis of symmetry.

CALCULATION OF THE OPTICAL FORCE AND THE TRAPPING POTENTIAL
The trap is characterized by the optical forces and the trapping potentials, which follow from the distribution of the electromagnetic (e.m.) field.We simulate this with Meep for positions on the axes shown in Figure 1 and using polystyrene beads as a force probe.The bead diameters are 0.2, 0.5, 1.0 and 1.4 µm.The optical force is evaluated by integrating the Maxwell stress tensor ↔ T over a cube with surface S around the bead: In Eq. ( 1) the brackets denote time averaging.Assuming summation convention, the Cartesian component F i (i = x, y, z) of the non-averaged force is written as Here n j is the j th component of the normal n to the surface S, while the tensor components T ij (i, j=x, y, z) are given by [12] Here E, D, H and B are the electrical field, the dielectric displacement field, the magnetic field and the magnetic induction, respectively; δ ij is the Kronecker delta and the asterix indicates complex conjugation.We use a cube of side a = 2 ( R + res), symmetrically placed around the bead, to evaluate the surface integral.Here R is the bead radius and res is the spatial resolution of the FDTD grid.For the e.m. fields at a certain time we use the values exported for two successive grid times, synchronized with Meep's synchronization tool [10].The time-averaged force results from averaging the timedependent force over the grid times in the last four wave periods of the simulation, where it has reached the steady state.We derive the trapping potential at position r 0 by integrating the time-average force according to To establish the grid size needed for sufficient accuracy, we perform simulation experiments.For the yand z-direction we take the calculated force for an empty cube (i.e.polystyrene replaced by water, so that the force should be zero) as a measure of the error in the force calculated for the particular bead size and resolution (in the experiments varied between 1 µm/10 and 1 µm/40).In these experiments the bead is positioned at a site of low symmetry, to avoid cancellation effects of integrals over symmetrically positioned cube faces.For the x-direction we perform experiments for res = 1 µm/20 and 1 µm/40, calculating the full force curve for the x-axis.In the trade-off between resolution and simulation time we arrive at res = 1 µm/20 for calculations on the x-axis and res = 1 µm/30 for calculations on the yand z-axis, respectively.For the x-direction this gives an underestimate of the amplitude of the force oscillations (see Section 5) of 35%, while in the yand z-direction this gives an estimated error in the force of about 20%.

FORCE CHARACTERISTICS OF THE DUAL BEAM TRAP
On the x-axis we see a composite behaviour: strong oscillations superimposed on a weakly sloping background line.The oscillations period is 297 nm, equal to the period in Figure 3(b), so that the oscillations arise from interference.We argue that in the gap a 1D optical lattice is formed, equivalent to an optical lattice formed with interfering laser beams and used to study trapped atoms [13].For bead diameters of 0.2, 0.5, 0.8 and 1.0 µm the background has the usual negative slope of the scattering force, while the slope increases with increasing bead size.At x = 0 the oscillatory curve for 0.2 and 0.8 µm has a negative slope, in agreement with a local gradient force pushing the particles to the intensity maximum.For 0.5 and 1.0 µm, on the contrary, the slope is positive and the particles are pushed to an intensity minimum.This behaviour, including the spatial periodic character, agrees with the results of Zemánek et al. [14], who theoretically study polystyrene beads in a 1D optical lattice in water.These authors find that the force depends periodically on the particle size, i.e. the particle's equilibrium position alternates between intensity maxima and minima, depending on size.As for the oscillations, the 1.4 µm bead falls in the same class as the 0.5 and 1.0 µm beads.The background, on the contrary, has a positive slope.
To analyse this in detail, we obtain force curves for bead sizes in the range 0.9-1.4µm, showing that with increasing size the negative slope weakens and turns into a positive slope.With increasing size, the bead increasingly probes the side of the intensity maximum close to the facet (Figure 2(c)).The µ maximum plays the role of a laser-beam focus and thus exerts an attractive gradient force on the particle.Further, upon approaching the facet, the scattering force on beads larger than the waveguide cross section (1 µm), shows a weaker increase than on smaller beads, an effect arising from the beam divergence.The total effect is that the gradient force exceeds the scattering force for the larger beads, reflected in the sign reversal.
For the yand z-direction the force is linear around the centre and has the usual restoring character (Hookian force).With increasing distance, the force shows an extremum, its position increasing with increasing bead size.The extremum occurs for a position such that the whole bead is on one side of the trap centre.This size dependent extremum indicates operation outside the Rayleigh regime (bead diameter << λ).
Dependence on bead size is consistent with the influence of the particle on the field seen in Figure 3.
The independence of the global intensity distribution on the relative phase of the beams (Section 3), implies that the global force distribution shows this property as well.The phase of the force oscillations for the x-direction, however, does depend on the relative phase.In a real device this enables moving the particle in the x-direction by changing the relative phase.Since the waveguides can be driven via a Y-junction (Section 7), a phase shifter in one of the waveguides will enable this particle motion.
The trapping force is strongest in the yand z-direction, where it takes maxima of 13.5 and 16 pN/W for the 1.4 µm bead.
From symmetry, one would expect these values to be the same.The difference indicates the accuracy of the calculations for res = 1 µm/30 (Section 4).Using these numbers, the average trap stiffness in the transverse is direction is 0.03 pN/nm/W, a value very suitable for trapping of µm-sized particles and comparable to that of the fiber-based dual-beam trap in [4].

TRAPPING POTENTIALS
To obtain the trapping potentials on the axes we integrate the force curves.Using Eq. ( 4), for the y-direction the trapping potential U y (y 0 ) at position y 0 (y 0 ≥ 0; the negative branch follows from symmetry) is given by A similar equation holds for U z (z 0 ).For the x-direction the force field is limited to the gap, while the size of the probe particle limits this range even further.In calculating U x (x) we integrate over the range for which the forces are displayed in Figure 4. We normalize U x (x) by applying a scaling factor which makes U x (0) equal to U z (0).
In Figure 5 we show the resulting trapping potentials U x (x) and U y (y), normalized to kT (T = 300 K).Since U y (y) is similar to U z (z), we omit it.For the 1.4 µm bead U x (x) is not shown, in view of the deviating slope of the background force, as discussed above.
U x (x), plotted in Figure 5(a) for P = 100 mW, shows a parabolic shape.The curvature increases with increasing bead size, while the oscillations correspond to the force oscillations.Stable optical trapping without bouncing at the channel walls can be induced for this power.This follows from the criterion U/kT >10 for stable trapping [15].The horizontal dashed lines in each potential correspond to U/kT = 10.For bead sizes 1.0, 0.8 and 0.5 µm the lines indicate fulfillment of the criterion.The average displacement of the trapped particle from the potential minimum is relatively large, as deduced from the crossings of the lines with the potential.For the 0.2 µm bead stable optical trapping does not occur for this power.Instead, for this size the channel walls take care of the confinement.At higher laser powers, where the depth of the minima of the oscillatory potential can exceed 10 kT, trapping may occur in a local minimum, similar to the situation of particle trapping in an optical lattice [13].
For U y (y), plotted in Figure 5(b) for P = 15 mW, the potential is deep enough to give stable trapping for 0.8, 1.0 and 1.4 µm beads, as again indicated by the dashed lines.Now stable trapping occurs for powers similar to those used in single beam laser traps.A difference with a laser trap is that our trap allows a larger average displacement from the potential minimum, which results from beam spreading instead of beam focussing.The power level of 15 mW for stable trapping is comparable to power levels used in Raman spectroscopy.Depending on the ratio of the power delivered to the particle to the input power, a higher input power may be needed to restrict the integration time of the spectrometer.It follows that the trap can also be used for Raman spectroscopy.

DISCUSSION AND SUMMARY
In a possible fabrication scheme, the waveguides are patterned in a Si 3 N 4 layer deposited on a glass wafer, yet omitting etching of the gap.The waveguides connect to a Y-junction patterned in the same layer and to be driven by a laser.For embedding, an SiO 2 overlayer is deposited and is then planarized.The channel and the gap are formed in a single patterning step.Bonding a second wafer on top, into which access holes to the channel are etched, seals the channel.Multiple Y-junctions enable multiple traps, by repeated branching out.
To have more design freedom to tailor the beams than with Si 3 N 4 cladded with SiO 2 , other platforms may be considered, e.g.SiON technology [16] and TriPleX technology [17].
In summary, we propose and simulate a novel dualwaveguide trap, for operation at near-visible and visible wavelengths.The trap is based on Si 3 N 4 waveguides with SiO 2 cladding.From our results it is expected that the device can trap biological objects such as bacteria.Fabrication of the trap is compatible with standard glass-based microfluidic technology.The maximum gradient force in the radial direction is 16 pN/W, corresponding to a trap stiffness of 0.03 pN/nm/W.This is comparable with values for fiberbased dual-beam traps, which cannot be produced with standard technology, contrary to our design.The depth of the trapping potentials indicates that our design yields stable trapping for practical optical powers.For coherent beams, we find an interesting oscillatory force in the axial direction, which can be used to axially move a particle.We envision that our design, when scaled to a lab-on-a-chip with multiple traps and equipped with an on-chip Raman spectrometer, may fulfill requirements for on-line sensing of drinking water and medical diagnostics.

FIG. 2
FIG. 2 Simulated intensity distribution of the electrical field in the xy-plane for the waveguide (a), and transverse (b) and axial (c) intensity profile for this waveguide, as further detailed in the text.In (a) the normalization is such that the full color scale is used.

FIG. 3
FIG.3Averaged intensity distribution for the complete dual-beam trap, in planes as indicated, with a 1 µm diameter polystyrene bead centered at point (0, 1, 0).The normalization of the distributions is such that the full color scale is used.

MFIG. 4
FIG. 4 Calculated force component as a function of position on the x-, y-, and z-axis, respectively (a-c), for several bead diameters.The insets in (b, c) are blow ups of the force curves for positive positions for the 0.2 and 0.5 µm diameter beads.The curves in (a-c) are piecewise cubic spline interpolations between the data points.

MFIG. 5
FIG. 5 Trapping potentials for the x-direction (a) and y-direction (b), as deduced from the curves in Figure 4.Note that the applied optical power is 100 and 15 mW in (a) and (b), respectively, while stable trapping can still be induced, as indicated by the position of the horizontal dashed lines.These are 10 kT above the bottom of the corresponding potential.