2-1. Theories
2-1-1. Optical force
In order to understand the behavior of particles around the waveguide, it is necessary to calculate the optical force acting on a particle. Assumptions and numerical analysis are used according to the particle size. For particles with diameter much larger than the wavelength of scattered light, the optical forces acting on a spherical particle were presented by Barton et al [14]. Mie’s scattering theory [15] and ray optics [16] have been employed as simple and excellent approximations that suffice for describing the behavior of the system.
As the particle have a radius ɑ which is much smaller than the wavelength of the scattered light (ɑ << λ/20), which is so called the Rayleigh particle, can be considered as with an induced dipole moment density P(r,t)
. Therefore, the electromagnetic force density can be described as the following equation:
P( ,t)
( ,t) ( ,t) ( ,t)12
The permittivity in the medium ism nm20. Using the vector identity and Maxwell equation, we can rewrite Eq. (2.1) as Poynting vector, the time-average modulus of the Poynting vector is the intensity I(r), which describes the number of photons per time and area in units of volt-amperes per square meter:
2
Photons carry this energy as well as a momentum m per volume at velocity c/nm . The momentum density (Abraham’s form) given in units of nanoseconds per cubic meter is
t)]
13
gradient force and a scattering force due to a strong field distribution near the focal point of a beam. However, these approximations fail to describe the forces when the particles and the waveguides are similar in dimension to the wavelength of the light used, and particles would experience significant distribution of evanescent fields. In this case, Maxwell stress tensor was a rigorous approach for calculating the optical forces.
The gradient and the scattering forces on a particle can be described as an electromagnetic stress using the Maxwell stress tensor:
TM DE HB (D E H B)I 2
1
(2.6)
where TM represents the Maxwell stress tensor, E is the electric field, B is the magnetic flux density, D is the electric displacement, H is the magnetic field, and I is the isotropic tensor.
Since we are interested in the transport processes which are much longer than the optical period on time scale, we use the time-averaged Maxwell stress tensor:
when expended out, the equation becomes
14
where the subscripts x, y, and z signify the coordinate directions. By integrating the time-averaged Maxwell stress tensor on the surface enclosing the particle, we can determine the total electromagnetic, FEM, force acting on the system by
FEM
S ( TM n)dS (2.9)2-1-1-2. Hydrodynamic analysis
Since the optical forces acting on particles are independent of the hydrodynamic properties of the fluid, the problems of optical and hydrodynamics can be decoupled. Here we use Navier-Stokes’s theory with continuity equations to describe the motion of particle in fluid under the influence of optical forces. At steady state, these equations are
v
is the flowing velocity of fluid, η is viscosity of the fluid, and P
is the pressure in
the system. Under these conditions the flow stress tensor TF
where I is the identity tensor, The hydrodynamic drag force FD on the particle resulting from fluid flow is described by
15
is surface normal vector.
Generally, hydrodynamic drag force is opposite to the optical force in direction.
Optical force balanced by fluidic drag force makes particle move with a constant speed, which is called terminal speed and is linearly proportional to the optical force. We use Navier-Stokes law to find the terminal velocity Vt for a particle of radius ɑ. For the simplest case, the particle is assumed to move through an infinite quiescent fluid. The hydrodynamic drag force FD is following modified expression.
t
D aCV
F 6
(2.15)
where C is a constant that contains the relationship between hydrodynamic drag force and the velocity. Fixen’s Law gives the constant term as,
where h is the distance between the center of the particle and the surface.
16
2-2. Finite element method
The finite element method (FEM) is a numerical technique for finding solutions of partial differential equations (PDE). For complicated optical systems, it can solve the boundary value problem, eigenvalue problem and find the steady state solution of a system by employing variational method. To apply this method, it requires discretizing a continuous domain into a set of discrete sub-domains, usually called elements, and the solution of each element would be approximated by certain characteristic form to solve the problems.
Here, the wave equations in the frequency domain for the magnetic and the electric fields are where c is the vacuum speed of light. To solve the wave equation for either magnetic or the electric field in frequency domain together with boundary conditions, standard FEM method proceeds in three steps. First, the wave equations are identified as solutions of certain variational problems where boundary conditions at the surface ∂V have been incorporated as additional terms of lagrangian L. The most general variational formulation for the electric field is given by
where μ is the magnetic permeability and ε is dielectric function both may varied in space. In addition, n
denotes the outward normal at the surface V and the electric field has to
17
satisfy the Dirichlet boundary condition n E
= 0 on ∂S. γe and U
are known quantities which are used to represent various other types of boundary conditions such as impedance boundary conditions and Sommerfeld radiation conditions. Finally, radiation sources within the computational domain V are described through the spatially varying current density j
.
Fig.2-1. The model (a) before (b) after mesh generation.
The second step is the most demanding step which consists of the discretization of the Lagrangian. The computational domain V is divided into a number of small-volume elements, the so-called finite elements. Within each element, the electric field is expanded into a series of certain elementary functions with unknown coefficients. It becomes possible to approximately enforce the div-conditions of the electric field within a given element as long as the dielectric function does not vary within this element. Fig.2-1 shows the model with/without mesh generation in someone case.
In the final step, these expansions facilitate the transformation of the variational eqs.
2.19 into a sparse set of linear equations via the Galerkin method. This matrix system can subsequently be solved via advanced linear algebra methods, either for obtaining eigenfrequencies and eigenmodes of the system of interest or to determine scattering cross sections of complex structures as well as transmittance and reflectance through functional elements.
18
2-3. Simulation results
Fig.2-2 illustrates the model system and computational domain used in this study. The model represents a silicon nitride waveguide which is on glass substrate and clad in water environment. As the evanescent mode from this hybrid system appeared, the suspend particles will be attracted by the extended field in the water medium, and thus the particle can access the mode and interact with the waveguide. Generally, the optical force can be divided into trapping and propulsion forces, which attracts and propels the particle, respectively.
Fig.2-2. The schematic of modeling system. A particle near the waveguide would experience a combination of trapping and propulsion force.
To characterize the trapping situation as well as the maximum trapping force in all three dimensions in our system, we used the commercial software package (COMSOL) to carry out the detailed three-dimensional finite element numerical analysis. The material properties of the waveguide structures, the glass substrate, and the surrounding water medium were taken into consideration while solving for the electromagnetic field distribution.
To determine the force exerted on a particle, a virtual spherical surface which enclosed
19
the entire particle was necessary to calculate the electromagnetic field on this surface. In this work, we assume the particle is separated from waveguide surface by 60nm and 1W input power is launched as the initial condition. The waves were launched into these waveguides with the boundary mode, in other words, we do not to take into account the contribution of the coupling loss. By using Maxwell stress tensor and integrating it over the surface enclosing particle, then we are able to obtain the optical forces exerting on a particle [17].
2-3-1. Waveguide geometry dependence of optical force
In order to investigating the dependence of optical force on waveguide geometry, we build the waveguide structures with different cross-sections here. Geometric parameters for the simulations used in our investigation are list in Table.2-1. The polystyrene particle with 1 and 2 μm in diameter dispersed in water can mimic most biological cells. We launch the transverse-electric (TE) polarized fundamental mode into each one to transport particles.
domain material Domain size Refractive index
Waveguide Silicon nitride width(w)×thickness(t) 2.2 (w = 0.5, 1 and 2μm,
t = 200 and 300nm)
Substrate Silicon oxide lower subdomain 1.45
cladding Water indefinite 1.33
particle polystyrene 1μm, 2μm 1.59
Table.2-1. Simulation and geometric parameters
Fig.2-3 illustrates wavelength dependence of the propulsion force (Fx) for a 2μm polystyrene sphere which is displaced at the horizontal center and separated from waveguide surface by 60nm. For each curve (waveguide of w = 0.5 μm and t = 200nm for example), the propulsion force first raises along with red-tuning the wavelength (region I) because the
20
evanescent field with longer wavelength extends deeper into water surrounding and interacts stronger with the particle. Then the force reaches a peak value when the interaction becomes the strongest. Furthermore, as the wavelength keeps red-tuning, the force drops gradually (region III), it is resulted from evanescent field extending beyond the particle and thus less evanescent field interacts with the particle. We use wavelengths λ of 1000nm, 1500nm and 2000nm to describe these profiles. Fig.2-4 shows the y component of electric field distribution at y = yparticle center with the corresponding wavelengths.
Fig.2-3. Wavelength-dependent propulsion force Fx for these geometries.
Fig.2-4. The distribution of y component of electric field at y=yparticle center for (a)λ=1000nm (b) λ=1500nm(c) λ=2000nm. Inset shows the streamline of power flow, respectively.
I
II III
21
When launching wavelength is 1000nm, the y component of electric field distribution is more condensed in the waveguide than others. Thus the particle experiences less interaction with the evanescent field. For the condition with 2000nm launching wavelength, although the amplitude of the y component of electric field is relatively small, the extend part of evanescent field interacts sufficiently with the particle. The propulsion force became larger than the previous condition. However, the losses also increased because of the more scattered wave. Here the most suitable wavelength is about 1550nm. Its evanescent field can provide the greatest interaction with the trapped particle.
Fig.2-5. Diagram showing the ray optics of a spherical Mie particle trapped and propelled on the waveguide in water surrounding. The evanescent wave driven forces is separated by propulsion force and downward trapping force.
These profiles can also be explained by the momentum transfer in the system. The proposed model is shown in Fig.2-5. The principal part of the momentum transfer from incident light to the particle is due to the scattered rays which are refracted by the particle.
The force exerted on a particle is given by the total change in the momentum of the incident photons. Because the momentum change of the ray from A to A’, we get a combination force which consisted trapping and propulsion forces. It also shows that a larger scattered angle θ
Fx
Fz
22
will result in stronger propulsion component of optical force on the particle in x direction.
Fig.2-5 also illustrates that the trapping force -Fz will be greater than the propulsion force Fx.
The streamline of the power flow have been showed in the inset of Fig.2-4. It is obviously that the incident light wave is scattered by the particle with the largest scattered angle when the λ=1500nm compared with other conditions, thus there is the strongest propulsion force.
Fig.2-3 also suggests that the thinner and narrower waveguide can reach stronger peak propulsion force at shorter wavelength. Base on this concept, we would be able to build an efficient transport system on a chip. For better coupling efficiency in experiments, we choose the thickness of 300nm instead of 200nm as our main parameter.
Fig.2-6. The diagram shows the wavelength-dependent transmissions for these waveguides in the single trapped particle system.
For the purpose of parallel manipulate, the guide wave must be strong enough to interact with other particles attracted subsequently after passing by one particle. We thus have to take into consideration that the transmission declines due to the existence of trapped
23
particle. We include this consideration by multiplying the propulsion force with transmission (T). The transmission is the percentage of guided power transmitted through waveguide section nearest to a single trapped particle. Fig.2-6 shows the transmission as functions of launching wavelength. It declines as red-tuning the launching wavelength.
The Fig.2-7 illustrates the propulsion force-transmission products as functions of launching wavelength. The high product guarantees that the waveguide will be able to transport a large amount of particles efficiently at the same time. Each peak wavelength is shorter than its counterpart in Fig.2-3. It means that we may sacrifice the propulsion force slightly for higher force-transmission product.
Fig.2-7. The products of propulsion force and transmission as functions of launching wavelength.
Similarly, we calculate the wavelength-dependent downward trapping forces Fz for these waveguides as shown in the Fig.2-8. The downward trapping force has a valley which represents the particle will be trapped most stably at the corresponding wavelength for each curve. It is worthy to note that, for each curve, the peak wavelength of the product almost coincides with the valley wavelength of the downward trapping forces. Hence we can
24
transport many particles simultaneously along the waveguide efficiently without easily disturbed when operated at the peak wavelength of force-transmission product.
Fig.2-8.Calculation of the downward trapping forces Fz on different width of waveguides for a wide range of wavelength with 300nm thickness.
Among these waveguides, 0.5μm wide waveguide is the most efficient and stable one.
As we operate the launching wavelength which near 1500nm, the better downward trapping and propulsion forces exerting on the particle are 614 and 142 pN, respectively. At the same time, the forces are much stronger than what wider waveguides provide near 1500nm in wavelength. Here we also calculate the resultant of forces of gravity and buoyancy is 2.05×10-3pN, and it is not at all affected for stable trapping.
25
2-3-2. Particle position dependence of optical forces
Fig.2-9. The diagram of cross-section of the system
To investigate stability of the trapped particle, we construct the system as shown in Fig.2-9. By the numerical calculation, we can detail the three-dimensional electromagnetic forces acting upon particle as functions of the position of particle at 1500nm and 1550nm in wavelength, which are upper and lower limits of the available wavelength in our measurement, respectively. From the force analysis, we were able to extract a parameter commonly used to describe optical traps: the trapping stiffness S which is defined as [18]
m equilibriu i
i
x S ( F )
(2-20)
it is the derivative of the restoring force with respect to the position of particle perturbed around the equilibrium point. We also note that the decay length (defined as the length where the force declines to e-1 times) for the force decreases in the z direction.
26
Fig.2-10. Comparison of the downward trapping forces Fz on a particle with 2μm in diameter using different waveguides and launching wavelengths (1500 and 1550 nm) at varying vertical positions.
Z=60nm λ=1500nm λ=1550nm
Width(μm) 0.5 1 2 0.5 1 2
Z-axis stiffness (pN nm-1W-1)
4.1 3.8 2.7 3.8 3.7 2.9
Decay length(nm)
143 105 88 153 111 92
Table.2-2. The characteristic trapping stiffness and decay length of the force profiles in Z-direction.
We calculate the forces at different positions in the lateral directions (perpendicular to x direction). These values are computed for 1W. We first compute the downward trapping force along z direction (Fz) at horizontal center of the waveguide as shown in Fig.2-10. A feature we observe is that the maximum Fz of these waveguides located at the position nearest to the waveguide surface. Fz declines fast once the particle separates from the waveguide surface for a few hundreds of nanometer. The trapping stiffness and decay length of the force profiles are shown in Table.2-2. The trapping stiffness for the 0.5μm wide waveguide is an
27
order of magnitude higher than that of silicon waveguide [19] (1.34pN nm-1 W-1 for 2000nm particle). And it is also larger than the stiffnesses for wider (1μm and 2μm) waveguides made of silicon nitride.
Fig.2-11. Comparison of the propulsion forces Fx at varying vertical positions.
Fig.2-12. The distribution of y component of electric field at y=yparticle center for (a) 0.5 μm (b) 1μm (c) 2μm wide waveguides.
We can also get the propulsion force distributed along z direction at horizontal center of the waveguide as shown in Fig.2-11. It similarly depends on width. The narrower waveguide have larger propulsion force. This situation can be explained in Fig.2-12. When the wavelength is operated at 1500 nm, the narrower waveguide have the stronger evanescent
28
field extending into water surrounding. The 0.5 μm wide waveguide exhibited the highest stiffness and propulsion forces.
Fig.2-13. Comparison of the transverse trapping forces Fy on a particle with 2μm in diameter using different waveguide at λ=1500 and λ=1550 nm.
Fig.2-13 illustrates the surface transverse trapping force (Fy) as a function of y position when particle is 60nm above the waveguide surface. Fy is weaker at the horizontal center of the waveguide, and becomes stronger when the particle is near the edge of the waveguide. For the same particle size, the narrower waveguide exhibits stronger Fy. This is due to the discontinuity of electric field at the edge of the waveguide, which is necessary for maintaining the continuity of electric displacement.
Y=0nm λ=1500nm λ=1550nm
Width(μm) 0.5 1 2 0.5 1 2
Y-axis Stiffness
(pN nm-1W-1) 0.34 0.21 0.08 0.35 0.22 0.09 Table.2-3. The characteristic trapping stiffness of the force profiles in Y-axis
The characteristic trapping stiffness of the force profiles are shown in Table.2-3. The
29
transverse trapping stiffness for the 0.5μm wide waveguide is higher than that of silicon waveguide (0.12pN nm-1 W-1 for 2000nm particle). And it is also larger than the stiffnesses for wider (1μm and 2μm) waveguides made of silicon nitride.
Fig.2-14. Comparison of the Propulsion forces Fx using different waveguide at λ=1500 and λ=1550 nm.
The dependence of propulsion force Fx on the position along y-axis is shown in Fig.2-14. It displayed the 0.5 μm wide waveguide can propel the particles with a wide range of transverse distance more efficient than the other waveguides.
30
2-3-3. Enhancement of optical forces by Indium Tin Oxide layer
In order to further improve the efficiency and the stability of particle manipulation, here we introduce an indium tin oxide layer between the waveguide and the glass substrate.
The schematic is shown in Fig.2-15.
Fig.2-15. The cross-section of the system with the ITO layer.
Because of the conductor-like optical property of ITO layer with complex dielectric function (εr= -2-0.5i near 1500nm in wavelength [20]), the guide mode will be forbidden extending into ITO film, and thus it distributes asymmetric in z direction with stronger evanescent field extending into water surrounding as shown in Fig.2-16(d) (with ITO film) compared to Fig.2-16(e) (without ITO film). A thicker ITO layer can force the guided mode to extend deeper into the upper surrounding as shown from Fig.2-16(d) to (a). The resulting optical force (Fx for example) upon 2 μm particle on 0.5 μm wide waveguide as a function of the thickness of ITO layer is shown in Fig.2-16.
31
Fig.2-16. Comparison of guided power distributions on 0.5 μm wide waveguide for (a) 300 nm (b) 200 nm (c) 100 nm (d) 30 nm (e)none thicknesses of ITO layer.
The indium tin oxide layer is conductor-like film, so the electromagnetic field is
The indium tin oxide layer is conductor-like film, so the electromagnetic field is