Simulations of Radiation PressureExperimentsbyMaximilien Bethune-WaddellB.A.Sc., The University of British Columbia, 2013A THESIS SUBMITTED IN PARTIAL FULFILLMENT OF THEREQUIREMENTS FOR THE DEGREE OFMASTER OF APPLIED SCIENCEinTHE COLLEGE OF GRADUATE STUDIES(Electrical Engineering)THE UNIVERSITY OF BRITISH COLUMBIA(Okanagan)August 2016c© Maximilien Bethune-Waddell, 2016i The undersigned certify that they have read, and recommend to the College of Graduate Studies for acceptance, a thesis entitled: Simulations of Radiation Pressure Experiments Submitted by Maximilien Bethune-Waddell in partial fulfillment of the requirements of The degree of MASc-Electrical Engineering . Dr. Kenneth Chau Supervisor, Professor (please print name and faculty/school above the line) Dr. Joshua Brinkerhoff Supervisory Committee Member, Professor (please print name and faculty/school in the line above) Dr. Jonathan Holzman Supervisory Committee Member, Professor (please print name and faculty/school in the line above) Dr. Ahmed Idris University Examiner, Professor (please print name and faculty/school in the line above) Dr. Murray Neuman External Examiner, Professor (please print name and university in the line above) August 10th 2016 (Date submitted to Grad Studies) Additional Committee Members include: Please print name and faculty/school in the line above Please print name and faculty/school in the line above AbstractElectromagnetic radiation, such as light, can exert a force on objects. This phenomenais referred to as radiation pressure. A popular example would be the radiation pressurefrom the sun which is presently the only force pushing NASA’s Kepler II satellite out of oursolar system. The most important material property relevant to radiation pressure is therefractive index n. The vacuum of space around the Kepler satellite has a refractive indexof n = 1. For this situation our ability to model radiation pressure is quite certain. Therehas been a debate over the last century on how to model radiation pressure in the presenceof matter, where the refractive index of the surrounding medium can be n > 1. Fiveprevalent models exist for electrodynamics, namely the Abraham, Minkowksi, Einstein-Laub, Chu, or the Amperian formulations. Various electrodynamic theories presented overthe last century have differences that only arise in exotic situations such as objects inliquid environments, material moving at relativistic speeds, and short timescales difficultto measure. Because technology has now advanced to potentially include these situations,this thesis addresses the renewed attention this debate requires. Each of these models aretested against two criteria using a simulation tool that solves for both electromagnetic andfluid dynamic phenomena. First, their compliance with the conservation laws of energy,momentum, and center-of-mass velocity. Second, their accord with experimental results. Asimulation environment is used to calculate the conserved quantities and observable effectsthat each model predicts under various experimental conditions. The simulation tool isa Matlab script that allows us to consider and compare the results of many experimentssimultaneously. Five significant experiments are analyzed in this thesis: the radiationpressure observed on metallic and dielectric mirrors, the deformation of a water-air interface,the deformation of a fluid-fluid interface, the displacement of polystyrene beads submergedin water, and the displacement of an oil droplet on a water surface. This thesis reachesthe conclusion that not enough data is available from past experiments to verify a singleelectrodynamic theory. Our work suggests simple new experiments that could.iiiPrefaceThis thesis was written under the supervision of Dr. Kenneth Chau, who co-authoredthe content and planned the analysis. The majority of this thesis represents the workdone to create a published article on radiation pressure authored by Dr.Chau and me [1].The waveguide experiments in Section 4.1 were carried out by Saimoom Ferdous and Dr.Thomas Johnson, a publication is in progress based on the results of their experiments. Twoparagraphs within Section 4.1 that describe and discuss the experiments were originallydrafted by Saimoom Ferdous and Dr. Thomas Johnson.ivTable of ContentsAbstract . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . iiiPreface . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . ivTable of Contents . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . viList of Figures . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . viiList of Tables . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . xiAcknowledgments . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 1Chapter 1: Introduction . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 1Chapter 2: Theory . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 52.1 Maxwell’s Equations . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 62.2 Energy Continuity and the Poynting Theorem . . . . . . . . . . . . . . . . . 72.3 Momentum Continuity . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 92.3.1 Example Calculation of Force Using the Stress Tensor . . . . . . . . 132.4 Center-of-Mass Velocity Conservation . . . . . . . . . . . . . . . . . . . . . 15Chapter 3: Methods . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 163.1 Finite-Difference Methods . . . . . . . . . . . . . . . . . . . . . . . . . . . . 173.1.1 The FDTD Method . . . . . . . . . . . . . . . . . . . . . . . . . . . 183.1.2 One-Dimensional Force Discretization and Field Averaging . . . . . 203.1.3 Two-Dimensional Force Discretization . . . . . . . . . . . . . . . . . 213.2 Finite-Difference Methods in Fluid Dynamics . . . . . . . . . . . . . . . . . 26Chapter 4: Results . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 274.1 The Balazs Thought Experiment . . . . . . . . . . . . . . . . . . . . . . . . 284.1.1 Momentum and Center-of-Mass Calculations . . . . . . . . . . . . . 304.1.2 An Experimental Variant on the Balazs Thought Experiment . . . . 344.2 Discussion and Conclusions on the Balazs Thought Experiment . . . . . . . 364.3 Momentum and Center-of-Mass Conservation For All Five Forms . . . . . . 37v4.4 Conservation in Two Dimensions . . . . . . . . . . . . . . . . . . . . . . . . 404.5 Radiation Pressure on Submerged Mirrors . . . . . . . . . . . . . . . . . . . 434.6 Radiation Pressure at Fluid Interfaces . . . . . . . . . . . . . . . . . . . . . 554.7 Optical Tweezing . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 624.8 Radiation Pressure in Negative Index Materials . . . . . . . . . . . . . . . . 684.8.1 Optical Tweezing With a Negative Index Bead . . . . . . . . . . . . 73Chapter 5: Conclusion . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 74Appendix A: Speeding up Matlab Code . . . . . . . . . . . . . . . . . . . . . . 77Bibliography . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 80viList of FiguresFigure 2.1 Demonstration of a pulse entering a lossy medium and leaving energybehind in the form of work done on the internal charges. . . . . . . . . . . . 8Figure 2.2 Differential volume element where the components of stress and shearare displayed for the plane defined with unit normal +xˆ. . . . . . . . . . . . 12Figure 2.3 A piece of steel immersed in a magnetic field ~B feels a force propor-tional to the spatial gradient of the magnetic field. . . . . . . . . . . . . . . 13Figure 2.4 Force versus separation z, normalized to the magnet length L. Valuesagree with simulated results except at very close proximity, where permanentmagnet field modeling requires more advanced methods [51]. . . . . . . . . . 14Figure 3.1 Figure showing how the derivative of f(x) can be expressed by thefinite difference between two points. . . . . . . . . . . . . . . . . . . . . . . 17Figure 3.2 The numerical error produced when field averaging is not implemented. 22Figure 3.3 Figure showing the discretized locations of the electromagnetic fieldand force values on the Yee grid. . . . . . . . . . . . . . . . . . . . . . . . . 22Figure 3.4 The general algorithm to solve for fluid deformation under radiationpressure. . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 27Figure 4.1 Illustration of the Balazs thought experiment. . . . . . . . . . . . . . 28Figure 4.2 A plane wave at the interface between free space and an impedance-matched medium (µr = εr). The plane wave is shown in terms of ~D and~B fields to illustrate the discontinuous nature of the Minkowski momentumdensity. The free space wavevector is described by k1 and propagation withinthe impedance-matched medium is described by k2. . . . . . . . . . . . . . . 30Figure 4.3 Energy density distributions at three moments in time predicted bythe Abraham (a) and the Minkowski (b) postulates, for the case of a pulsenormally incident onto a slab under conditions identical to those in Figure 4.4. 31Figure 4.4 Pulse momentum exchange with an impedance-matched medium. . . 33viiFigure 4.5 (a) Infrared image at 4 minutes showing the temperature at points (1)and (2) on either end of the Teflon segment.(b) The waveguide experimentalsetup where a 13.56 MHz source propagates a TEM wave in the ~k directionalong an impedance-matched waveguide that is terminated with a matchedload. . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 35Figure 4.6 Salient dimensions for the impedance-matched waveguide labeled onthe waveguide with their values tabulated. This configuration representsa close analogue to the impedance-matched medium of the Balazs thoughtexperiment. . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 36Figure 4.7 Simulations of the Minkowski, Abraham, Einstein-Laub, Amperian,and Chu postulates with conservation of momentum and center-of-mass ve-locity. . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 38Figure 4.8 Instantaneous energy density distributions at three moments in timepredicted by (a) the Abraham, Einstein-Laub, and Chu postulates, (b) theMinkowski postulates, and (c) the Amperian postulates, for the case of apulse normally incident onto a slab under conditions identical to those inFigure 4.7. . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 39Figure 4.9 Two-dimensional field-based implementation of the Balazs thoughtexperiment. . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 41Figure 4.10 Non-degenerate center-of-mass displacement predictions made by theMinkowski (MN), Abraham (AB), Einstein-Laub (EL), Chu, and Amperian(AMP) postulates for closed systems in which finite pulses are obliquelyincident onto tilted impedance-matched slabs having (a) n =√3 (blue) and(b) n˜ = −1.3 + i0.07 (yellow). . . . . . . . . . . . . . . . . . . . . . . . . . . 42Figure 4.11 Momentum transfer between a pulse (λ0 = 632 nm, τ = 2 fs) in anenclosure containing a fluid (n = 1.6) and a submerged mirror approximatedas a perfect electric conductor. . . . . . . . . . . . . . . . . . . . . . . . . . 45Figure 4.12 Momentum transfer between a pulse (λ0 = 632 nm, τ = 2 fs) in anenclosure containing a fluid (n = 1.6) and a submerged mirror approximatedas a perfect magnetic conductor. . . . . . . . . . . . . . . . . . . . . . . . . 47Figure 4.13 Key assumptions used to model the Jones-Richards experiment. . . . 49Figure 4.14 Normalized time-averaged radiation pressure on submerged mirrorspredicted by the Minkowski, Abraham, Einstein-Laub, Chu, and Amperianpostulates. . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 51viiiFigure 4.15 Time-averaged force density distributions predicted by the Minkowski,Abraham, Einstein-Laub, Chu, and Amperian postulates, shown for the caseof a continuous-wave beam (λ0 = 632.8 nm) normally incident onto a di-electric mirror composed of alternating λ0/4-thick layers of MgF2 and ZnSsubmerged in a fluid of refractive index n = 1.6. . . . . . . . . . . . . . . . 53Figure 4.16 Normalized radiation pressure on a submerged dielectric mirror pre-dicted by the Einstein-Laub, Chu, and Amperian postulates as a functionof the force density integration distance into the fluid for cases in which thedielectric mirror is capped by either the lower-index layer (blue solid) or thehigher-index layer (red dashed). . . . . . . . . . . . . . . . . . . . . . . . . . 54Figure 4.17 Instantaneous (top row) and time-averaged (bottom row) electro-magnetic force density distribution exerted by a continuous-wave beam (λ0= 632.8 nm) normally incident from air onto water (n = 1.33) predicted by(a) the Abraham (solid blue) and Minkowski (dashed red) postulates and(b) the Einstein-Laub and Chu postulates (solid black) and the Amperianpostulates (dashed orange). . . . . . . . . . . . . . . . . . . . . . . . . . . . 57Figure 4.18 Two-dimensional fluid dynamic simulations representative of the Ashkin-Dziedzic experiment [15]. The lower fluid (blue) represents water (n=1.33)and the upper fluid (white) represents air (n=1). . . . . . . . . . . . . . . 59Figure 4.19 Two-dimensional fluid dynamic simulations representative of the ex-periments by Casner et al. [19, 20]. . . . . . . . . . . . . . . . . . . . . . . . 61Figure 4.20 Simulated radiation pressure on a polystyrene bead (green) placed atthe edge of the focal point of a continuous-wave beam (λ0 = 532 nm) in abackground medium of water (n = 1.33). . . . . . . . . . . . . . . . . . . . 64Figure 4.21 Simulated dynamics of an air bubble placed at the edge of the focalpoint of a TM-polarized, continuous-wave beam (λ0 = 532 nm) in a back-ground medium of water (n = 1.33, σ = 0.7 mN/m, η = 1 mPa·s). . . . . . 65Figure 4.22 Tractor beaming of an oil droplet (n = 1.42) resting at the interfacebetween air (n = 1.0) and water (n = 1.33). . . . . . . . . . . . . . . . . . . 66Figure 4.23 Tractor beaming of a dielectric particle using two obliquely incidentplane waves. . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 67Figure 4.24 A perfectly absorbing enclosure is partially filled with a negative-indexslab. . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 69Figure 4.25 A perfectly absorbing enclosure is filled with negative index fluid andan Ag mirror. . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 71ixFigure 4.26 Demonstration of the momentum a mirror experiences with respectto its surrounding refractive index. All forms predict degenerate results fora refractive index n > 1, but differences arise for n < 1. . . . . . . . . . . . 72Figure 4.27 Optical tweezing of an impedance matched negative index bead underTM illumination. . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 73Figure 5.1 A minimum working example (MWE) to demonstrate how to usehandles to update a figure as opposed to re-calling the “figure” function.The code above will produce a figure of a simple sine wave. . . . . . . . . . 78Figure 5.2 An example to demonstrate how vectorizing for loops improves thespeed of matlab code. . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 79xList of TablesTable 1.1 : Momentum density, stress tensor, and power flux corresponding to the Minkowski,Abraham, Einstein-Laub, Amperian, and Chu formulations of electrodynamics.xiAcknowledgmentsForemost, I would like to acknowledge Dr. Kenneth Chau for his unwaivering support,motivational character, dedication to quality, and for being a good example as a humanbeing. The mindset and skills I have developed by working under Ken have proven my mostvaluable assets and will surely continue to shape my professional and personal life for thebetter. Ken and I had a lot of fun during this work exploring different aspects of scienceranging from the fundamentals of quantum mechanics, the nature of absorption, the myster-ies of radiation pressure, and in developing numerical methods for electrodynamic modeling.As well, I would like to thank Dr. Joshua Brinkerhoff for his contribution to the fluiddynamics portion of this thesis and his role as a committee member. Also, a thank youto Dr. Jonathan Holzman for his thoughtful comments in reviewing this thesis and lettingme join his lab’s lunch meetings, where I have developed many friendships and gained ex-posure to teamwork driven experiments. I would additionally like to thank Dr. StephenO’Leary for his role in creating a Mitacs internship where I have applied the simulationtools of this thesis in an industry setting. Another acknowledgment to Ken, Holzman,and Stephen for their roles as academic supervisors during the various research projectsI have completed over the last five years during my undergraduate and graduate degrees.These projects have led to me being awarded an NSERC PGS-D scholarship. I will usethis scholarship to further my understanding of nature and contribute to reducing society’sentropy under a subsystem assumption. Lastly, I have enjoyed the good company of myfellow graduate students and want acknowledge the influence of Waqas Maqsood, ReyadMehfuz, Samuel Schaefer, Naomi Fredeen, Mohammed Al Shakhs, Asif Al Noor, BrandonBorn, Iman Aghanejad, and Chris Collier.1Chapter 1: IntroductionIn the early 1900s, works such as James Maxwell’s treatise [2] developed most of thetheory behind our present understanding of electromagnetic phenomena. One of the mostimportant results within Maxwell’s work was the proof that light was an electromagneticwave carrying momentum, meaning they could exert a force via Newton’s first law. Thisled Maxwell to predict light could produce a pressure on objects, a phenomenon now re-ferred to as radiation pressure. Shortly after Maxwell’s prediction, measurements of mirrordeflection due to incident radiation from a vacuum [3, 4, 5] confirmed the existence ofradiation pressure. In vacuum, mass-energy equivalence (E = mc2) and the definition oflinear momentum can be combined to show that the momentum of light is related to itsenergy content via p = E/c. Experiments in vacuum have confirmed that a mirror gainsmomentum of 2E/c due to the change in momentum upon reflection. In the presence ofmatter, competing theories have been introduced to model radiation pressure. Authors suchas Minkowski [6] and Abraham [7] began a now century-old debate over how to model themomentum density of light in matter. Minkowski proposed a form of momentum densitygiven by ~D× ~B, with ~D representing the displacement field and ~B the magnetic flux density.On the other hand, Abraham introduced a form given by ( ~E × ~H)/c2, with ~E represent-ing the electric field, ~H the magnetic field, and c the speed of light. Although identicalin vacuum, the Abraham and Minkowski forms give divergent predictions of the electro-magnetic momentum density in matter. If we consider an electromagnetic wave enteringa dielectric medium of refractive index n from free space, the Minkowski form predicts anelectromagnetic momentum increase proportional to n, whereas the Abraham form pre-dicts an electromagnetic momentum decrease proportional to 1/n. Dichotomy between thetheories of Abraham or Minkowski became known as the Abraham-Minkowski controversy.This debate has grown in scope over the past century to include alternative electrodynamictheories such as the Einstein-Laub, Amperian, and Chu formulations [8, 9, 10, 11, 12].2Among these electrodynamic theories, compelling claims have been made in supportof either the Minkowski or Abraham forms of electromagnetic momentum. Experimentalmeasurements of radiation pressure on submerged mirrors [13, 14], fluid interfaces [15, 16,17, 18, 19, 20], semiconductors carriers [21, 22, 23], and dilute atom gases [24] have consis-tently revealed a refractive-index proportionality most commonly interpreted to support theMinkowski momentum density [25, 26, 27, 28]. On the other hand, convincing theoreticalarguments have been developed in support of the Abraham momentum density. Balazs [29]put forward an extension of the Einstein photon-in-box thought experiment [30] to showthat only the Abraham momentum density preserves linear center-of-mass translation in aclosed system consisting of a rigid, reflection-less dielectric block and an electromagneticpulse. Center-of-mass translation, in combination with energy and momentum, is a quan-tity that must be conserved for a closed system [30, 31]. Strong evidence for the Minkowskiand Abraham momentum densities has led to an acceptance of the potential correctness ofboth [26, 32, 27, 33] and a heuristic approach to their application: The Minkowski or “wave”momentum can be used to predict radiation pressure phenomena, whereas the Abraham or“kinetic” momentum can be used to predict center-of-mass translations [34, 35, 36, 37, 38,39].In the 1960s, Penfield and Haus [11] showed that multiple electrodynamic formulationscan be consistent with global energy and momentum conservation by separating a closedsystem into sub-systems. For example, a closed system containing matter and electro-magnetic fields could be broken down into a kinetic sub-system describing the mechanicalmotion of matter, a thermodynamic sub-system describing the energy flow in matter, andan electrodynamic sub-system describing the momentum and energy stored in the fields.Each sub-system is modeled using energy and momentum continuity equations, which takeon the generic forms∇ · ~S + ∂∂tW = −φ, (1.1)and∇ · T¯ + ∂∂t~G = −~f. (1.2)where T¯ is the stress tensor, ~G the electromagnetic momentum density, ~f the force density,~S the energy flux density, W the energy density, and φ the power density. In a closedsystem, the continuity equations of the various sub-systems, indexed by the variable l, arecoupled so thatΣφl = 0, (1.3)3andΣ~fl = 0, (1.4)which means that energy and momentum move between systems via φl and ~fl to enablethe entire model to conserve energy and momentum.Attempts were made in the mid-1900’s [13, 14, 15] to test different postulates of momen-tum by measuring radiation pressure on illuminated objects, but consistent interpretationof these experiments has been elusive. Radiation pressure on submerged mirrors [14], forexample, has been explained in terms of the Minkowski momentum [25, 26, 27, 28], theAbraham momentum [40], the arithmetic mean of the two [41, 42, 43, 44], or a canonicalmomentum of quantum mechanical origins [45]. To reconcile the disparate literature, twomutually-exclusive hypotheses have emerged: either the definition of momentum is totallyarbitrary and there are many simultaneously valid electrodynamic theories [33, 46, 47, 39]or there is a single theory that uniquely predicts observable phenomena but has escapedexperimental detection [48]. This thesis presents evidence to argue the latter.A consolidated analysis is performed on the five electrodynamics theories mentioned sofar (namely: Abraham (AB), Minkowski (MN), Amperian (AMP), Einstein-Laub (EL), andChu) to investigate their theoretical and experimental consistency. Table 1.1 shows eachof their respective combinations of stress tensor, power flux density, and momentum density.4FormMomentumDensityPower Flux Stress TensorMinkowski† ~D × ~B ~E × ~H ( ~D · ~E + ~B · ~H)I¯/2− ~D ~E − ~B ~HAbraham ~E × ~H/c2 ~E × ~H [( ~D · ~E + ~B · ~H)I¯ − ~D ~E − ~B ~H − ~E ~D − ~H ~B]/2Einstein-Laub~E × ~H/c2 ~E × ~H (εo ~E2 + µo ~H2)I¯/2− ~D ~E − ~B ~HAmperian εo ~E × ~B ~E × ~B/µo (εo ~E2 + µ−1o ~B2)I¯/2− εo ~E ~E − µo−1 ~B ~BChu ~E × ~H/c2 ~E × ~H (εo ~E2 + µo ~H2)I¯/2− εo ~E ~E − µo ~H ~HTable 1.1 . Momentum density, stress tensor, and power flux corresponding tothe Minkowski, Abraham, Einstein-Laub, Amperian, and Chu formulations ofelectrodynamics. †The Minkowski energy-momentum tensor is presented in itsclassical asymmetric form [6, 49] as well as it’s symmetric form using a powerflux of c2 ~D × ~B [8].In this thesis, Section outlines the fundamental theory of how energy and momentumconservation is modeled for electromagnetic field quantities. Section uses numerical meth-ods to solve for field and force density, and explains them in detail, as they are used topresent the results. These methods are implemented using finite-difference techniques. Us-ing this solver, results are presented as numerical calculation of various important radiationpressure experiments conducted to date. First, we re-visit Balazs’ thought experiment inSection 4.1 to demonstrate its bias towards favoring the Abraham momentum density us-ing arguments of center-of-mass translation. Following this, numerical methods are furtherused in Sections 4.5- 4.8 to recreate the results from experiments on fluid-submerged mir-rors, radiation-induced deformation of fluid interfaces, the manipulation of small dielectricparticles by laser beams, and a thought experiment concerning radiation pressure in neg-ative index media. Our calculations show that all five postulates can be made consistentwith conservation laws; however, not all forms are consistent with experimental results. Inparticular, the postulates proposed by Abraham and Einstein-Laub are shown to have themost consistency with all radiation pressure experiments considered.5Chapter 2: Theory2.1 Maxwell’s EquationsAs we are concerned with electromagnetic fields, the most fundamental equations toconsider are Maxwell’s equations, namely: Ampere’s law, Faraday’s law, Gauss’s law, andGauss’s law for magnetics, given by 2.1 to 2.4.∇× ~H = ∂∂t~D + ~J, (2.1)∇× ~E = − ∂∂t~B + ~Jm1, (2.2)∇ · ~D = ρ, (2.3)and,∇ · ~B = 0, (2.4)respectively, where ~E is the electric field, ~H the magnetic field strength, ~J the currentdensity, and ~Jm the magnetic current density. These equations completely describe thespatial and temporal evolution of electromagnetic field quantities in the presence of sources.The force on a test charge q can be described by the Lorentz force densityf = q ~E + q~v × ~B, (2.5)where the charge q experiences force components due to electric and magnetic fields. Be-cause the force contribution of the magnetic field component is perpendicular to the direc-tion of particle motion, the magnetic field can not do work on moving charges. It will beshown in later sections that this is different than situations where magnetic fields can dowork such as in the case of magnetic materials. The description of motion resulting fromelectromagnetic influences is the study of electrodynamics. A correct electrodynamic theory1does not exist, but provides a mathmatical tool for useful features numerically, for example a perfectlymatched layer.6must obey conservation of energy, momentum, and center-of-mass velocity. Conservationof energy for an electromagnetic system is mathematically stated as∫v(∇ · ~S)dv +∫v(∂∂tW )dv = −∫Uwdv, (2.6)where ~S is the power flux density in [W/m2], Wtot is the energy density, and Uw is the workdone between the fields and any free charges. Any electrodynamic theory must also obeyconservation of momentum∫v(∇ · T¯ )dv +∫v(∂∂t~G)dv = −∫~fdv, (2.7)where T¯ is the stress tensor [Pa], ~G the electromagnetic momentum density [Kg/(s ·m2)],and ~f the force density [N/m3]. All quantities in equations 2.6 and 2.7 are assumed hereto only be associated with interaction between fields and charges, meaning they will onlyconsist of field and current density values. However it should be noted these equationsare generalizable to many different systems. For example, the Navier-Stokes equation thatdescribes fluid dynamics is derived from momentum continuity very similar to equation 2.7.2.2 Energy Continuity and the Poynting TheoremEquation 2.6 is called Poynting’s theorem [50]. Derivation of ~S and W for an electro-magnetic system will typically start with an expression for the work done by the fields. Aspreviously mentioned only electric fields can do work on charges. Using this principle wecan generalize the work done, Uw, as a function of only the electric field ~E and the currentdensity ~J asUw =∫v( ~J · ~E)dv, (2.8)which can be inserted into equation 2.6 by integrating over an area Ac that encloses v togive∫Ac(~S · nˆ)dA+∫v(∂∂tW )dv = −∫v( ~J · ~E)dv. (2.9)We can further re-arrange equation 2.9 using Maxwell’s equations to arrive at expressionsfor W and ~S based only on field quantities. Substituting Ampere’s law and setting ~J =∇× ~H − ∂∂t ~D along with the vector identity ∇ · ( ~E× ~H) = ~H · (∇× ~E)− ~E · (∇× ~H) gives−∫v( ~J · ~E)dv =∫v(∇ · ( ~E × ~H) + ∂∂t12[ ~E · ~D + ~H · ~B])dv. (2.10)7Identification of the terms on the right side leads us to the conventional expression for thePoynting vector~S = ~E × ~H, (2.11)and energy density,W =12[ ~E · ~D + ~H · ~B]. (2.12)Equation 2.10 can now be used to describe wave dynamics in lossy systems. For example,Figure 2.1 depicts a pulse initially entering a system with 1 J and leaving behind half ofits energy due to damping associated with induced charge movement. This loss can becalculated either directly using equation 2.8, or with only field values using equation 2.12.The best way will depend on what numerical or analytical technique is being used. Formost numerical methods the fields ( ~E, ~D, ~B, and ~H) are more readily available than thecurrent density, making equation 2.12 a preferred choice to solve for the energy density.eee e1.0 J 0.5 J 0.5 J0 J 0.5 J 0.5 JFigure 2.1: Demonstration of a pulse entering a lossy medium and leaving energy behindin the form of work done on the internal charges.82.3 Momentum ContinuityWe next explore the implications of momentum continuity by examining equation 2.7∫Ac(∇ · T¯ )dA+∫v(∂∂t~Gtot)dv = −∫~fdv.The stress tensor, T¯ , is a second rank tensor and its divergence results in a vector withxˆ, yˆ, and zˆ components. The approach to derive the conventional forms of the Maxwellstress tensor T¯ is to start by equating the force density in equation 2.7 to the Lorentz forcedensity in equation 2.5,∫v~fdv =∫v(ρ ~E + ~J × ~B)dv. (2.13)We assume our expression for ~f is solely responsible for any mechanical momentum and canthen parse ~Gtot into mechanical and electromagnetic parts, ~Gtot = ~GEM + ~Gmech. For nowwe assume there is no mechanical component to T¯ but note that many models incorporatea view where T¯ = T¯EM + T¯mech and argue that the use of T¯mech is how different formsof momentum density can be proven to be equivalent [33]. For now we will assume onlyan electromagnetic component to the stress tensor. To derive T¯ and GEM one can plugρ = εo∇ · ~E and ~J = 1µo (∇× ~B)− εo ∂∂t ~E into equation 2.13 to get∫Ac(∇ · T¯ )dA+∫v(∂∂t~GEM )dv = −∫(ddt~Gmech)dv = −∫v~fdv (2.14)and∫v~fdv =∫vεo[ ~E(∇· ~E)− c2( ~B×(∇× ~B)) + ~B × ∂∂t~E]dv.Next, we use the identity ~B × ∂∂t ~E = − ∂∂t ~E × ~B + ~E × ∂∂t ~B and add c2(∇· ~B) = 0 to get∫v~fdv =∫vεo[ ~E(∇· ~E)− ~E×(∇× ~E)− c2( ~B×(∇× ~B)) + c2 ~B(∇· ~B)− ∂∂t( ~E × ~B)]dv,which can be re-written as−∫vεo[ ~E(∇· ~E)− ~E×(∇× ~E) + c2 ~B(∇· ~~B)− c2( ~B×(∇× ~B))]dv + εoµo∫v∂∂t( ~E × ~H)dv= −∫v~fdv.(2.15)9From the left-hand side of equation 2.3, we can identify the momentum density~GEM =1c2( ~E × ~H)and in order to view our expression as a continuity equation we need to re-phrase the otherterms to fit into a stress tensor∇·T¯ = −∫vεo[ ~E(∇· ~E)− ~E×(∇× ~E) + c2 ~B(∇· ~B)− c2( ~B×(∇× ~B))]dv. (2.16)We consider only the x-component of the electric terms, knowing that similar analysis can becarried out for the other components. Development of equation 2.3 gives the xˆ-componentof[ ~E(∇· ~E)− ~E×(∇× ~E)]xˆ= Ex(∂∂xEx +∂∂yEy +∂∂zEz)xˆ− Ey( ∂∂xEy − ∂∂yEx) + Ez(∂∂zEx +∂∂xEz)xˆ(2.17)and now we add 0 = 12∂∂xE2x − 12 ∂∂xE2x and group the product rule terms to get∂∂xE2x +∂∂yExEy +∂∂zExEz − 12∂∂x(E2x + E2y + E2z ).This same form repeats for each Cartesian component to give[ ~E(∇· ~E)− ~E×(∇× ~E)]αˆ =∑β=[x,y,z]∂∂β(EαEβ − 12( ~E · ~E)δαβ) (2.18)andδαβ =1 if α = β0 else . (2.19)Equations 2.18 and 2.19 describe the divergence of a second rank tensor. When expanded,the notation for each component is given by∇·Tαβ = ∇·Txx Txy TxzTyx Tyy TyzTzx Tzy Tzz, (2.20)where, for example, the x components Txx, Txy, and Txz would be stated as10Txx = εo[E2x −12(E2x + E2y + E2z )]Txy = εo[ExEy]Txz = εo[ExEz]. (2.21)The total electric-field-dependent components for all terms in equation 2.3 are[E2x − 12 (E2x + E2y + E2z )] [ExEy] [ExEz][EyEx] [E2y − 12 (E2x + E2y + E2z )] [EyEz][EzEx] [EzEy] [E2z − 12 (E2x + E2y + E2z )] . (2.22)Through a similar derivation, the magnetic-field-dependent components of equation 2.3 are[B2x − 12 (B2x +B2y +B2z )] [BxBy] [BxBz][ByBx] [B2y − 12 (B2x +B2y +B2z )] [ByBz][BzBx] [BzBy] [B2z − 12 (B2x +B2y +B2z )] . (2.23)Each of these terms are the electromagnetic contributions to the components of shear, τ ,and normal force, σ, within a volume density. For example, Figure 2.2 shows a smalldifferential volume element with components of normal force, σij , and shear, τij . If weconsider a differential cube such as in Figure 2.2, each face of that particle has a normalforce and a shearing force. The first subscript, i, defines the normal axis plane, and thesecond subscript j refers to the force direction on that plane. Shown in Figure 2.2 are thecomponents associated with the plane normal to the x-axis, where the shear on this planein the zˆ direction is given by τxz, the shear in yˆ direction is given by τxy and the forcenormal to this plane is σxx. This treatment can be applied to all six faces of the particle.11σxxτxyτxzxˆyˆzˆ σxx τxy τxzτyx σyy τyzτzx τzy σzzFigure 2.2: Differential volume element where the components of stress and shear aredisplayed for the plane defined with unit normal +xˆ.122.3.1 Example Calculation of Force Using the Stress TensorEquation 2.7 is not used often in electromagnetics as a means to calculate force. Typ-ically, the Lorentz force equation 2.5 is used directly. To understand this equation, wecan consider it within a simple situation where electric and magnetic fields are not timevarying. In Figure 2.3, a piece of steel is immersed in a static magnetic field ~B. We wish tocalculate the attraction between a bar magnet and a piece of steel. There is no electric fieldpresent, to simplify analysis further we can assume the magnetic field is ~B = Bozˆ in freespace and Bo/µr zˆ in the steel which has relative permeability µr. Under these assumptions,equation 2.7 has no time derivative terms and becomes∫v(∇ · T¯ )dv = −∫v~fdv. (2.24)Force exists only where (∇ · T¯ ) is non-zero. This occurs at the interface between steeland air. Using the Bz dependent terms of equation 2.23, we can solve for the force as anintegral at the surface of the steel. If we consider the bottom face of the steel to havesurface area S′, the force on this surface is∫v(∇ · T¯ )dv =∫v−~fdv =∫v(− ∂∂zB2z2µo)dv = F zˆ.Assuming dz straddles the interface equally from −dz/2 to dz/2, the net force F isF =12µo[Bo2 − Bo2µr]S′. (2.25)steeldzBozˆS’Figure 2.3: A piece of steel immersed in a magnetic field ~B feels a force proportional to thespatial gradient of the magnetic field.13Equation 2.3.1 can be interpreted as a “gradient” force due to the spatial variation in~B. As well this force is associated with the work done in moving the steel a small dzwithin the external magnetic field [51]. This example is comparatively simpler then usingthe Lorentz force density, which would involve modeling the steel as a volume density ofcirculating current loops. Using equation 2.3.1, the force on steel from a permanent magnetas a function of distance is shown below in Figure 2.4. Experimental measurements of themagnet force data was obtained from a permanent magnet manufacturer [52].0 0.5 1 1.5050100z/LF[N]SimulatedK & J Magnetics Inc. [52]Figure 2.4: Force versus separation z, normalized to the magnet length L. Values agree withsimulated results except at very close proximity, where permanent magnet field modelingrequires more advanced methods [51].142.4 Center-of-Mass Velocity ConservationThe final fundamental requirement we shall consider for a good electrodynamic theoryis center-of-mass conservation. First, consider a system of two particles m1 and m2, movingat respective speeds of v1 and v2 in the xˆ direction. Conservation of linear momentum forthis system gives(m1~v1 +m2~v2) = Pxˆ, (2.26)where ~P is the total linear momentum of the two particles. The global center-of-massequation with respective particle locations of x1 and x2 would bem1x1 +m2x2 = (m1 +m2)x¯. (2.27)where x¯ is the system center-of-mass location. The center-of-mass velocity is then the timederivative of this equation,m1dx1dt+m2dx2dt= (m1 +m2)dx¯dt= P. (2.28)We have now made the left-hand side of equation 2.28 equal to the constant momentum ofthe system ~P , meaning that conservation of linear momentum will also yield conservationof the system center-of-mass velocity. This demonstrates that conservation of momentumand center-of-mass velocity are linked for a system of particles, as equation 2.26 was derivedfrom equation 2.28. However, these two conservation laws cannot be considered mutual forelectromagnetics because the definition of electromagnetic mass. For example the electro-magnetic mass of a pulse could be calculated asme =∫vW (x, y, z)/c2dv, (2.29)where W (x, y, z) is the spatial distribution of the electromagnetic energy density. Thisleaves the possibility that a spatial distribution of W will conserve global linear momentumwith total momentum me but not center-of-mass velocity. This idea has been consideredby authors such as Balazs in [29, 53]. As this is a discussion on mass and energy, Einstein’smass-energy equivalence E = mc2 can be stated for electromagnetic fields as~S = ~Gc2. (2.30)15This equivalence of power flux density ~S and momentum density ~G can be found in [54,pg. 27 §9], and in more detail related to fields in [31, pg. 260 §75]. Modern texts inelectrodynamics such as [55, pg. 285 sect. 11.11] and [56, pg. 269 sect 12.5] present thatpower flux density and momentum density are in general related by a factor of c2 in bothvacuum and matter. The implementation of equation 2.4 for a given ~G is said to be a“symmetric” consideration [57, 49, 31]. Although the relation of symmetry to center-of-mass conservation was derived by Einstein in relation to relativity [30], we will show inSection 4.1 that it affects predicted observables for classical phenomena.16Chapter 3: Methods3.1 Finite-Difference MethodsMany engineering problems are too complex to be done analytically. Computers haveenabled us to simulate these problems using numerical methods. Heat transfer, fluid me-chanics, and electromagnetics, each have a governing set of differential equations that canbe numerically solved using finite difference methods. The underlying principle behind thefinite difference method is shown in Figure 3.1. The derivative of a function f(x) can be ex-pressed as the “finite difference” between two points. A smaller step size, ∆x, means a moreaccurate finite difference. Shown is a plot of f(x) = x2, if an initial value of the function atsome xc is known, its neighboring values can be found via f(xc±∆x) ≈ f(xc)± f ′(xc)∆x.This is the basic concept behind all finite difference methods. The more formal definition ofa finite difference, and its accuracy, is related to the Taylor series expansion [58, pg. 104].0 0.5 1 1.5 2 2.5 3 3.50510f ′(2) ≈ ∆f∆x ≈ f(3)−f(1)3−1xf(x)1Figure 3.1: Figure showing how the derivative of f(x) can be expressed by the finite dif-ference between two points. As seen, the larger the step size between points, the moreinaccurate the finite difference.173.1.1 The FDTD MethodThe Finite-Difference-Time-Domain (FDTD) method is one of the most well developedfinite difference methods in computational electromagnetics. A full and complete accountof its implementation can be found in [59]. As the name suggests, it is a time-domaintechnique, the equations used are the time-dependent Ampere’s law 2.1 and Faraday’slaw 2.2. We will derive the FDTD equations in 1D, assuming a plane wave with only Hyand Ex components propagating in the zˆ direction. Amperes law simplifies to∂∂zHy = εo∂∂tEx.Next, we discretize space and time choosing index i for space and index n for time. Fornotation purposes we will assume spatial separation is expressed by ∆z and time separationby ∆t. We need to place the fields somewhere within these indexes, we will put Hy in timeat n+ 12 and take its spatial derivative at i,[∂∂zHn+ 12y]i=1∆z[Hn+ 12y (i+12)−Hn+ 12y (i− 12)]= εr(i)εo[∂∂tEx(i)]n+ 12.To be consistent with our notation taking the time derivative at i forces us to place thecentered difference of Hy between (i +12 ) and (i − 12 ). If assume Hy exists at n + 12 , theEx field must be located at n so that its centered difference is equal about n+12 . The fulldiscretized equation becomes1∆z[Hn+ 12y (i+12)−Hn+ 12y (i− 12)]=εr(i)εo∆t[En+1x (i)− Enx (i)],which can be re-arranged to solve for En+1x (i) as a function of past values. The same canbe done using Faraday’s law to re-arrange for Hn+ 12y (i+12 ). The notation of using half stepssuch as i + 12 is referred to as a “staggered grid” and is accredited to Yee [60]. The finalupdate equations areEn+1x (i) = En(i) +∆tεoεn+ 12r (i)∆z[Hn+ 12y (i+12)−Hn+ 12y (i− 12)](3.1)andHn+ 12y (i+12) = Hn−12 (i+12)− ∆tµoµnr (i+12 )∆z[Enx (i+ 1)− Enx (i)] , (3.2)where it has been intentionally emphasized that the material parameters of relative per-mittivity εn+ 12r (i) and permeability µnr (i+12 ) are also staggered in time and space. This isusually inconsequential as permittivity and permeability rarely vary in time and their stag-18gered spatial locations cause negligible errors. If they did vary in time, or a problem withvery high spatial variation was considered, more complicated methods than FDTD may haveto be employed, such as the Method of Moments, or adaptive mesh refinement (AMR) [61,62]. For stability, equations 3.1.1 and 3.1.1 require that the spatial and temporal step sizesbe related by∆t ≤ ∆z/(2c),where c is the speed of light. If ∆t were larger than this value, then the algorithm couldnot capture wave propagation over ∆z. The finite difference scheme would then be unableto approximate a solution to a wave equation.193.1.2 One-Dimensional Force Discretization and Field AveragingThis section will introduce how to calculate ~S, W , T¯ , and ~f using a one-dimensional(1D) example. If we consider plane wave propagation along the xˆ axis, the only componentsof the magnetic and electric fields are Hy and Ex. We can place them numerically at Enx (i)and Hn+ 12y (i +12 ), where i is the spatial index, and n is the temporal index. There is nostandard for placement of ~S, W , T¯ , and ~f numerically on the staggered Yee grid in termsof i and n. However, there is advantageous placement for reduced numerical error and fastcomputation. Averaging numerical quantities is quite common in other finite differencesolvers such as in fluid dynamics and heat transfer [63]. The force fz, via equation 2.7 is∂∂tGz +∂∂zTzz = −fz.We can then start our derivation with placing force numerically at fnz (i)1, express the otherterms as central differences, and re-arrange for fnz (i) asfnz (i) = −1∆z[Tnzz(i+12)− Tnzz(i−12)]−∆t[Gn+ 12z (i)−Gn−12z (i)]. (3.3)Where the act of placing ~f at the index location of i and n in equation 3.3 required thestress tensor to be located at Tnzz(i +12 ), and the momentum density at Gn+ 12z (i). As anexample, for the Abraham form these quantities are explicitly given byTnzz(i+12) = Dx(i)nEx(i)nBy(i)nHy(i)nand,Gn+ 12z (i) = En+ 12x (i)Hn+ 12y (i).Note the fields are not defined at By(i)n or Hy(i)n in the normal magnetic field updateequation 3.1.1. They can be solved by averaging from known quantities. For example,By(i)n can be approximated asBy(i)n =14[By(i+12)n+12 +By(i+12)n−12 +By(i− 12)n+12 +By(i− 12)n−12](3.4)and En+ 12x (i) as1We could also solve for the mechanical momentum, ~g directly through ∂∂tg = f , or numerically as[~gn+12 (i) − ~gn− 12 ]/∆t = ~fnz (i).20En+ 12x (i) =12[En+1x (i) + En−1x (i)]. (3.5)The discretization for W would beWn(i) = Wn−1(i)− ∆t∆z[Sn+ 12z (i+12)− Sn+ 12z (i− 12)], (3.6)which requires additional field averages. For example, Sn+ 12z (i+12 ) is given bySn+ 12z (i+12) = En+ 12x (i+12)Hn+ 12y (i+12)where En+ 12x (i+12 ) is a field quantity that does not exist in the normal FDTD Yee algorithm.It can be approximated by field averaging asEn+ 12x (i+12) =14[Ex(i+ 1)n+1 + Ex(i+ 1)n−1 + Ex(i− 1)n+1 + Ex(i− 1)n−1]It would be desired to reduce the number of times averaging has to be done, however thestaggered mesh algorithm of FDTD makes field averaging unavoidable.The benefits of field averaging is demonstrated in Figure 3.2 where an energy densitycalculation from an FDTD simulation is used as an example. In the figure a pulse entersthe simulation space and then passes through an impedance-matched slab. The slab islossless so a correct solution should show no energy density left behind. The figure showsthat using FDTD to track the energy density of a pulse without field averaging can causesignificant error, in the case shown, up to 2.5 %. When field averaging is used, the error isnegligible at the same step size.3.1.3 Two-Dimensional Force DiscretizationWe can expand our previous discussion on force discretization in one dimension to dis-cuss problems that require two dimensional (2D) analysis. To do this, we will considerhow to numerically solve for the Minkowski force density using the corresponding stresstensor, T¯MN and momentum density, ~D × ~B. We will consider the 2D staggered Yee gridshown below in Figure 3.3 where only Ex, Ey, and Hz field components are present. Thiscorresponds to “transverse-magnetic” (TM) propagation.210100W/Wo%02.5Field AverageingNo Field AverageingW/Wo%(a)(b) (c)t=18 fs t=35 fsFigure 3.2: The numerical error produced when field averaging is not implemented. A pulseis directed at an impedance-matched slab from free space with refractive index n = 2. (a)A DC error value of energy density is at both the entrance and exit faces of the slab afterthe pulse propagates through. In (b) this error becomes negligible when field averagingis implemented. (c) The pulse as it exits the slab with no field averaging, the error iscomparable in value to the field energy itself when no field averaging is done.Eny(i, j) (i+ 1, j)(i, j + 1) (i+ 1, j + 1)EnxHn+ 12zfn+ 12x( ij+12)(i, j) (i+ 1, j)(i, j + 1) (i+ 1, j + 1)fn+ 12y(i+12j )Wn+ 12(i+12j+12)Figure 3.3: Figure showing the discretized locations of the electromagnetic field and forcevalues on the Yee grid. The placement of fx and fy are such to reduce the number of fieldaverages that need to be taken, this reduces the numerical error.The total momentum ~G in continuity equation 2.7 can be parsed into electromagnetic andmechanical contributions where ~G = ~GMN + ~Gmech. The force density is then calculatedas22∇ · T¯MN + ∂∂t~GMN = − ∂∂t~Gmch = −~fMN , (3.7)where the Minkowski stress tensor is T¯MN =12 (~D · ~E + ~B · ~H)I¯ − ~D ~E − ~B ~H, and theresulting Minkowski force density is denoted as ~fMN . We now need to make a choice ofwhere to place each term on the Yee grid. The example presented here will only considerthe xˆ-component of force. The xˆ component of the momentum density term ∂∂t~Gx will beplaced spatially at(ij+ 12)and temporally at (n+ 12 ) to give∂∂tGn+ 12x( ij+12)=1∆t[(DyHz)n+1( ij+12)− (DyHz)n( ij+12)].We can see from the Yee grid in Figure 3.3 that Dy is already located correctly, but Hzis not spatially or temporally at the right location, to remedy this we can approximateHnz( ij+12)from the field averageHnz( ij+12)=14Hn+ 12z(i+12j+12)+Hn− 12z(i+12j+12)+Hn+ 12z(i− 12j+12)+Hn− 12z(i− 12j+12) . (3.8)By using averaging such as in equation 3.8 we can use much larger step sizes before numericalerror becomes apparent. To calculate the x component for force, the discretized form ofequation 3.7 can now be expressed asTn+ 12x( ij+12)+∂∂tGn+ 12x( ij+12)= −fn+ 12x( ij+12). (3.9)Which requires another field average to solve for Tn+ 12x( ij+12). Using the placement in equa-tion 3.1.3 we can express this placement asTn+ 12x( ij+12)=12[Tn+1x( ij+12)+ Tnx( ij+12)].The full expression for Tn+ 12x( ij+12)needs to be solved by numerically expanding(∇·T¯MN ) · xˆ = [ ∂∂xTxx +∂∂yTxy +∂∂zTxz]xˆ.For our TM example, these tensor components simplify toTxx =12[DxEx +DyEy +BzHz]−DxEx,Txy = −DxEy,23andTxz = 0.These expressions can be expressed as four different unique terms for discretization that wecan denote here by t1, t2, t3, and t4 within the expressionTnx( ij+12)=12[−tn1( ij+12)+ tn2( ij+12)+ tn3( ij+12)]− tn4( ij+12). (3.10)The terms in equation 3.10 expand totn1( ij+12)=1∆x(DxEx)n(i+12j+12)− (DxEx)n(i− 12j+12) ,tn2( ij+12)=12∆x[(DyEy)n( i+1j+12)− (DyEy)n( i−1j+12)],tn3( ij+12)=1∆x(BzHz)n(i+12j+12)− (BzHz)n(i− 12j+12) ,andtn4( ij+12)=1∆y[(DxEy)n( ij+1)− (DxEy)n(ij)].The field averages required for t1 areDnx(i+12j+12)=12[Dnx(i+12j )+Dnx(i+12j+1)]andDnx(i− 12j+12)=12[Dnx(i−12j )+Dnx(i−12j+1)].For t3 they areDnx(ij)=12[Dnx(i+12j )+Dnx(i−12j )]andDnx( ij+1)=12[Dnx(i+12j+1)+Dnx(i−12j+1)].For t4 they areEny(ij)=12[Eny( ij− 12)+ Eny( ij+12)]andEny( ij+1)=12[Eny( ij+32)+ Eny( ij+12)].24The t2 term does not require any averaging as the Ey and Dy components already exist atthe (i ± 1, j ± 12 ) locations on the Yee mesh. This averaging process can be repeated foreach of the different momentum density forms and stress tensors listed in Table 1.1 . TheFDTD method provides a direct way to solve for electrodynamic quantities, but care mustbe taken in regards to field placement and field averaging.253.2 Finite-Difference Methods in Fluid DynamicsFor the fluid studies in this thesis, we implement a CFD method based on the markerand cell (MAC) front tracking technique, as outlined by Tryggvason [64]. A complete step-by-step guide meant for students can be found in [63]. The equation to be solved is theincompressible Navier-Stokes equationρ∂u∂tf= ρ∇ · uu−∇P + qo∇2u+ ρg + fo + 〈fe〉, (3.11)where qo∇2u is the diffusion term with viscosity qo, ρ∇ · uu the advection term, ∇P isthe pressure gradient, ρg is the effect of gravity, fo are body forces including surface ten-sion, and 〈fe〉 is the time-averaged electromagnetic force. The electromagnetic force is timeaveraged because the time scales for fluid deformation is much larger than the electromag-netic timescales considered in this thesis. Interaction between electromagnetic fields anddeformable fluids is described in a feedback loop cyclically modeling each on its relevanttime scales under quasi-static conditions. The time-averaged electromagnetic force densityis calculated using the FDTD method and acts as a body force within the MAC method,which then determines the deformation of the fluid under the combined influence of gravity,radiation pressure, surface tension, and viscous forces. Once the front interface separatingthe two fluids has displaced one coarse grid spacing ∆, the fluid boundaries are updated onthe Yee grid, and the time-averaged electromagnetic force density (averaged over at least 30cycles) is re-computed. In this way, the FDTD and MAC solvers progress in turn, feedingeach other updated distributions of the force density and fluid-interface location at the endof their turns. The result is a series of detailed snapshots of fluid deformation under varyingradiation pressure and material distribution. The algorithm for this continuous operationis shown in Figure 3.4. The initial conditions for the fluid velocity at all locations is zero.26fn+ 12y(i+12j )fn+ 12x( ij+12)Wn(i+12j+12)vq(i+12j )uq( ij+12)ρq+ 12(i+12j+12)fq′t(r)fn+ 12y(i+12j )fn+ 12x( ij+12)Wn(i+12j+12)vq(i+12j )uq( ij+12)ρq+ 12(i+12j+12)fq′t(r)fn+ 12y(i+12j )fn+ 12x( ij+12)Wn(i+12j+12)vq(i+12j )uq( ij+12)ρq+ 12(i+12j+12)fq′t(r)1(a)(b)Figure 3.4: (a) The general algorithm to solve for fluid deformation under radiation pres-sure. The density, ρ, and the material parameters of permittivity and permeability areinitialized together, with the electromagnetic grid having a finer grid resolution that thenfluid grid. The electromagnetic fields are then solved and the resulting force density ~fe istime-averaged before inserted into the Navier-Stokes algorithm. (b) staggered grid locationsfor the electromagnetic values of fx, fy, and energy density W along with staggered gridlocation for the fluid dynamic values of the velocity components u, v, and density ρ. Alsoshown are the force contributions from the fluid surface tension ft, which exists on a frontwith sub-cell resolution relative to both staggered grids.27Chapter 4: Results4.1 The Balazs Thought ExperimentThe Balazs thought experiment was a theoretical argument made in 1953 that comparedthe Abraham and Minkowski momentum densities on the premise of center-of-mass veloc-ity conservation [29]. In his experiment, Balazs compares two hypothetical closed systems.Both systems contain a transparent impedance matched slab of length L with mass M anda pulse of mass m << M . In the first system the pulse flies over the slab and no interac-tion occurs, the total system momentum is mc and the system’s center-of-mass velocity ismc/(m + M). In the second system the pulse collides with the slab and slows down to aspeed of c/n. During this interaction, the total momentum of each system must be equalto mc as there are no external sources of work to differentiate them [29] the pulse withinthe interacting system will have lagged behind the first system’s pulse. Figure 4.1 depictsboth systems of the Balazs thought experiment with lag d.M(a)mpota tb(b)xˆLdFigure 4.1: Illustration of the Balazs thought experiment. (a) In the first system a pulse ofmass m does not interact with the slab, meaning global momentum is mc and global centerof mass velocity is mc/(M +m). (b) In the second system where the pulse enters the slaband slows down to c/n.28We will now re-produce the derivation of this thought experiment similar to the authorsin [65, 66, 67] and show the expanded steps used by Barnett [66]. We take the electro-magnetic energy of a pulse to be that of a photon h¯ω, and its mass given by m = h¯ω/c2.With the slab initially stationary, the global conserved momentum of the system at ta isequal to only the photon momentum mc = h¯ω/c. Taking the time interval in Figure 4.1 as∆t = tb − ta, the difference in the photon’s center of mass location between system 1 and2 is(h¯ωc∆t)−(h¯ωcn∆t)=h¯ωc∆t(1− 1n),where the time interval can be expressed as ∆t = nLc to giveh¯ωc2L (n− 1) . (4.1)For uniform center-of-mass conservation between the two systems, the slab contribution tocenter-of-mass within the interacting system in Figure 4.1 (b) must be equal to equation 4.1.Equating the product of M∆xs to equation 4.1 yields a slab displacement of∆xs =h¯ωMc2L (n− 1) . (4.2)This displacement can be used in the equation for momentum conservation in the secondsystem and leads to an expression for the photon momentum within the slab. Momentumconservation for the second system can be expressed ash¯ωc= M∆xs∆t+ pem, (4.3)where equation 4.2 represents ∆xs which can be combined with ∆t =Lnc to givepem =h¯ωcn, (4.4)with equation 4.4 corresponding to the Abraham momentum density. For the remainderof this thesis we will consider center-of-mass and momentum calculations in a field baseddescription where the electrodynamic interactions are solved over continuous time intervals.This will provide useful insight beyond considering a particle picture and comparing eventsonly in terms of before and after an interaction. For example we will be able to discusshow the Minkowski form could still be considered valid in terms of center-of-mass despitethe derivation just shown.294.1.1 Momentum and Center-of-Mass CalculationsAs opposed to the assumption of a fixed pulse mass m, we can combine equations 2.6and 2.4 to express an electromagnetic mass density in terms of any postulated ~G. We willconsider a case where there is no external work (Uw = 0), giving an energy continuityequation of∇ · ~S + ∂∂tW = 0,where expressing the electromagnetic mass as m = W/c2, and substituting ~S = ~Gc2 givesm(t) =∫ t0−(∇ · ~G)dt, (4.5)where ~G can be any form of momentum density. The Minkowski momentum density, ~D× ~B,is made of discontinuous fields, which means that implementing the Minkowski momentumdensity yields a surface energy current that can be solved for by taking the volume integralof equation 1.1 at the interface between two different media. This is analogous to the currentdensity that results when deriving the boundary conditions for the electric displacement field~D [68]. We can demonstrate this with a simple 1D-plane-wave example, as in Figure 4.2,where the interface between free space and a lossless impedance-matched medium (εr = µr)is considered.zˆxˆyˆεo, µo~D1 = εoaosin(k1z − ωt)xˆ~B1 =µoaoηosin(k1z − ωt)yˆεr, µr~D2 = εoεraosin(k2z − ωt)xˆ~B2 =µoµraoηosin(k2z − ωt)yˆδ/2l/21Figure 4.2: A plane wave at the interface between free space and an impedance-matchedmedium (µr = εr). The plane wave is shown in terms of ~D and ~B fields to illustrate thediscontinuous nature of the Minkowski momentum density. The free space wavevector isdescribed by k1 and propagation within the impedance-matched medium is described byk2.30When implementing the energy continuity equation 2.6 at the boundary in Figure 4.2with δ → 0, a surface energy current Kw(t) can be defined as~Kw(t) =∫ t0(~S1 − ~S2)dt, (4.6)where ~S1 is the power flux density in free space and ~S2 represents the power flux densityinside the medium. This expression can be non-zero as the fields associated with theMinkowski power density can be discontinuous, as shown in Figure 4.2. The magnitude ofthe surface energy density [J/m2] associated with the case in Figure 4.2 isKw(t) =a2o2ηo[(1− εrµr)t+ (εrµr − 1)(sin2(2ωt)2ω)]. (4.7)For a case more analogous to the pulse in the Balazs thought experiment, Figure 4.3illustrates the electromagnetic energy density distribution as a pulse interacts with animpedance matched lossless slab of index (n = 2). Values are calculated using the FDTDmethod in Section . For the symmetric Minkowski form in Figure 4.3, illumination of theslab causes an energy extraction at the entrance face and s deposition of energy at the exitface.AB MNenergydensity(a)(b)z (µm)t=17 fs t=28 fs t=39 fsFigure 4.3: Energy density distributions at three moments in time predicted by the Abra-ham (a) and the Minkowski (b) postulates, for the case of a pulse normally incident ontoa slab under conditions identical to those in Figure 4.4. The discontinuous nature of theMinkowski energy density predicts surface energy spikes analogous to those calculated inFigure 4.2.31For the same slab as in Figure 4.3, Figure 4.4 compares momentum and center-of-masscalculations for the Abraham and Minkowski postulates. As shown in the snap-shots of theelectric field, the pulse is created at the left end of the enclosure and propagates to the righttowards the slab. During this interaction, global momentum and center-of-mass velocity areconserved (i.e., both remain fixed at zero). After the pulse enters the slab, its reduction inforward momentum is compensated by the speeding up of the slab and the increase in slabmomentum. However, if we replace the Abraham momentum density with the Minkowskimomentum density while leaving everything else unchanged, center-of-mass conservation islost (dashed line in Figure 4.4 (d)). Upon pulse entry, the increase in the forward pulsemomentum is balanced out by the backward slab momentum, which causes the systemcenter-of-mass velocity to stray from zero. Violation of this conservation principle has beenone of the primary arguments against the validity of the Minkowski momentum density [29,27, 37]. Using the symmetric form of the Minkowski power flux (solid line) allows the pulsemass to change in spatial density upon entering the slab and conserve the global center-of-mass velocity. This pulse mass affects the slab mass, as energy continuity implies that anychanges in the pulse mass must be balanced by local changes in the slab mass.32−Po0Po/nPoAB−Po02nPo MN-∆z00∆z0derpalerpABpulse: ~E × ~Htotal: ~E × ~H10 15 20 25 30 35 40 45-∆z00∆z0time (fs)derpalerpMNt=6 fs t=17 fs t=28 fs t=39 fsnormalizedmomentumnormalizedz¯pulse slabenclosure total(a)(b)(c)(d)Figure 4.4: Pulse momentum exchange with an impedance-matched medium. The resultsare normalized to the initial free space momentum of Po, and slab index n. (a) The Abrahammodel predicts a pulse momentum decrease proportional to 1/n when immersed in the slab,where as the Minkowski model (b) predicts a pulse momentum increase proportional to n.(c) Conservation of center-of-mass velocity for the Abraham form. (d) Conservation ofcenter-of-mass velocity for the Minkowski form. The Balazs implementation for the pulseand total center-of-mass displacement is shown with a dotted line, using the Abrahampower flux density of ~E× ~H, the Minkowski form does not conserve center-of-mass velocity.The solid line in (d) illustrates that when paired with its symmetrical form of power fluxdensity, c2 ~D × ~B, the Minkowski form can conserve center-of-mass velocity.334.1.2 An Experimental Variant on the Balazs Thought ExperimentHypothesizing the existence of a surface energy such as in equation 4.6, we can considerwhat would happen if this energy density were real and able to thermalize. To model this,we can relate equation 4.6 to a volume generation term in the one-dimensional heat diffusionequation,∂2T∂z2+q˙(z, t)k=1α∂T∂t, (4.8)where q˙(z, t) is the internal energy generation, k is the thermal conductivity, and α is ther-mal diffusivity.We propose the following experiment to distinguish the Abraham and Minkowski formsof momentum density. At microwave frequencies, a coaxial guided wave structure propagat-ing a dominant TEM mode can easily be constructed [69]. Within the coaxial structure, adielectric material is inserted in a section of the guide to create a boundary condition wherepower flows through a uniform medium. This concept is illustrated in Figure 4.6. Sincethe characteristic impedance of the transmission line is determined by the conductor radii,a and b, and the dielectric constant of the medium, the inner radius can be stepped downto maintain a constant characteristic impedance throughout the length of the transmissionline [69, ch. 5]. If all sections of the transmission line are matched, then no reflections atthe dielectric boundaries are created and a TEM wave should propagate without appre-ciable loss. The thermal properties for the Teflon used in this experiment are consideredat their room temperature values. Teflon thermal conductivity k is 0.25 W/(m ·K), andα is 0.124 mm2/s [70]. According to Figure 4.3, a cooling effect may be expected at theentrance and a heating effect would exist at the exit. If a portion of the surface energyfor the Minkowski form is assumed to thermalize, a lower bound on the amount of energyconverted from electromagnetic to thermal should be proportional to the loss tangent ofthe material. We set the total surface generation density to be equal to δ ·Kw(t) calculatednumerically using equation 4.6 where δ is the loss tangent of Teflon, reported by [71] tobe 1.2× 10−6 at 13.56 MHz. This low loss tangent preferentially creates a situation whereno observed heating would indicate a Abraham form, but an observation of appreciableheating suggests the Minkowski form.34A picture of the coaxial structure is shown in Figure 4.6, and consists of a Telfondielectric in the mid section sandwiched between two air dielectric regions. The coaxial lineis designed to have a characteristic impedance of 50 Ω and measurements are made in ahigh-power test bed where the incident source power can be adjusted to several hundredwatts. A frequency of 13.56 MHz was used with an off-the-shelf 1 kW power amplifier asa source to create large incident power fluxes. The coaxial line is terminated in a high-power 500 W load which is matched to 50 Ω. Teflon has a dielectric constant of 2.1 atthe operating frequency of 13.56 MHz; therefore, very low insertion loss is expected. Twosmall holes, labeled as spots 1 and 2 in Figure 4.5 are drilled into the outer conductor tocreate apertures where temperature can be imaged with an infrared camera. Figure 4.5 (b)shows a temperature difference of 0.4 K, which is much smaller than the predicted 2.0 Ktemperature difference from solving equation 4.8. Most importantly, Figure 4.3 predicts theexit face to be hotter than the entrance face relative to the propagation direction, a trait thatholds for continuous wave illumination as well. This contradicts the observed temperaturedifference where Spot 2 in Figure 4.5 (a) is actually cooler than the entrance face of Spot1. We conclude that the absence of measurable heating or cooling at the interfaces of theblock more directly supports the Abraham momentum and Poynting vector without relyingon whether the arguments of Laue are correct [49].Spot 1 Spot 2Transmission Line(b)(a) Matched Load~kFigure 4.5: (a) Infrared image at 4 minutes showing the temperature at points (1) and(2) on either end of the Teflon segment.(b) The waveguide experimental setup where a13.56 MHz source propagates a TEM wave in the ~k direction along an impedance-matchedwaveguide that is terminated with a matched load.35Air εo, µoL2L1Teflon εr, µr2a1 2bCu2a21~ka1 = 1.59 mma2 = 1.29 mmb1 = 4.345 mmL1 = 151.7 mmL2 = 69.88 mmFigure 4.6: Salient dimensions for the impedance-matched waveguide labeled on the waveg-uide with their values tabulated. This configuration represents a close analogue to theimpedance-matched medium of the Balazs thought experiment.4.2 Discussion and Conclusions on the Balazs ThoughtExperimentIs the surface energy density predicted by the Minkowski postulates real and physical?The possibility of similar discontinuities in energy density have been discussed at the sur-face of magnetizable media when invoking an Amperian definition of power flux ~E × ~Bµo.It has been argued that these energy discontinuities are “hidden” and purely mathematicaland not physically measurable [72]. Hidden energy in the Amperian definition of powerflux has been rationalized because the fundamental building blocks of magnetization arecurrent loops, which possess inherent energy necessary to sustain perpetual current motion.Similar hidden energy arguments have yet to be put forward to reconcile the Minkowskipower flux, although these arguments would likely invoke hidden energy storage in bothcurrent loops and polarization charge.We have presented here an alternative route to comparing the Abraham and Minkowskimomentum forms using an experimental implementation of the Balazs thought experiment.We have examined their differences in predicted energy deposition at the interface of losslessdielectrics. The results of our experiments show that the symmetric form of the Minkowskitensor produces a surface energy term. In light of the fact that the Abraham form hasshown no requirement for the existence of an anomalous surface energy current to conserveall variables, we agree with Balazs’ conclusion that the Abraham form is more likely to becorrect.364.3 Momentum and Center-of-Mass Conservation ForAll Five FormsThe Balazs thought experiment compared only the Minkowski and Abraham forms ofmomentum density; however, we can also use it to study the Amperian, Einstein-Laub, andChu formulations. Shown below in Figure 4.7 is the momentum and center-of-mass conser-vation for all five electrodynamic theories for the more general case of a dispersive slab withloss. Here, the Einstein-Laub and Chu formulations are identical to the Abraham form,which is to be expected as the Einstein-Laub and Chu stress tensors in Table 1.1 are onlydifferent within their two-dimensional shearing terms. The Amperian form exhibits someunique properties. The Amperian form predicts a decrease in electromagnetic momentumcompensated by an increase in slab momentum, having the same directional behavior asthe Abraham form, but different in magnitude. The center-of-mass conservation for theAmperian form demonstrates (dotted line) that pairing the Amperian momentum densitywith the canonical Poynting vector ~E× ~H also causes a violation of center-of-mass velocityanalogous to how the Minkowski form violates center-of-mass conservation shown in Fig-ure 4.4. Therefore, also analogously when the power flux of the Amperian form is solvedfor by using equation 2.4, it does conserve center-of-mass velocity.370p0Pulse SlabEnclosure Total0∆z0derpalerp0p00∆z0Pulse: ~E × ~HTotal: ~E × ~H0p00∆z0Pulse: ~E × ~HTotal: ~E × ~H0 50 0 50time (fs) time (fs)t=18 fst=0 fs t=26 fs t=37 fsnormalizedcenterofmassdisplacementnormalizedmomentum(a)(b)(c)(d)(e)(f)AB, EL, & ChuMNAMPAB, EL, & ChuMNAMPFigure 4.7: Simulations of the Minkowski, Abraham, Einstein-Laub, Amperian, and Chupostulates with conservation of momentum and center-of-mass velocity. Shown here isa general case of a lossy dispersive impedance-matched slab with refractive index n¯ =2.04 + i0.01, using electric and magnetic drude parameters of γe = γm = 2.5 × 1014,ωpe = ωpm = 1.2 × 1015, and ∞ = µ∞ = 2.2. The pulse is centered at a free-spacewavelength λ0 = 632 nm with a pulse width τ = 2 fs. Both the pulse and slab are immersedin vacuum and contained in a rigid enclosure. (top row) Snapshots of the pulse electric fieldemitting from the left side of the enclosure, propagating through the slab, and disappearingon the right side of the enclosure. (a)-(c) Momentum contained in the pulse (red), slab(blue), and enclosure (black) as described by (a) the Abraham, Einstein-Laub, and Chupostulates, (b) the Minkowski postulates, and (c) the Amperian postulates. Momentum isnormalized to the free-space momentum of the pulse, p0. The total momentum (green) isalways conserved and fixed at zero. (d)-(f) Center-of-mass displacement of the pulse (red),slab (blue), and enclosure (black) as described by (d) the Abraham, Einstein-Laub, andChu postulates, (e) the Minkowski postulates using either ~E× ~H (dashed) or ~D× ~B (solid)as the power flux, and (f) the Amperian postulates using either ~E× ~H (dashed) or 0 ~E× ~B(solid) as the power flux. The center-of-mass displacement is normalized to the net slabdisplacement divided by its mass ratio. The system center-of-mass (solid-green) predictedby the Minkowski and Amperian postulates with the Poynting vector ~E× ~H (dashed green)strays from zero and thus, is inconsistent with center-of-mass velocity conservation.38The energy density for the momentum and center-of-mass calculations in Figure 4.7 isshown in Figure 4.8. As with the previous analysis in Figure 4.3, surface energy densityexists for the discontinuous Minkowski and Amperian power flux forms. In addition, someenergy is left behind in the bulk of the material for all forms proportional to the loss of thedispersive slab.energydensityz (µm)AB, EL, & Chu MN AMP(a)(b)(c)Figure 4.8: Instantaneous energy density distributions at three moments in time predictedby (a) the Abraham, Einstein-Laub, and Chu postulates, (b) the Minkowski postulates, and(c) the Amperian postulates, for the case of a pulse normally incident onto a slab underconditions identical to those in Figure 4.7.394.4 Conservation in Two DimensionsWe next model the electrodynamic interaction between a pulse and a slab in two-dimensions. Extension of the analysis to more than one dimension is important, as itincorporates additional sources of electromagnetic momentum associated with electrostric-tion, magnetostriction, field gradients, and out-of-plane scattering (refraction) which areunavoidable in many experiments. As well, the shear terms previously mentioned for thestress tensor of the Einstein-Laub and Chu formulations become applicable. We use all fivesets of electrodynamic postulates to study two configurations in which a finite-width pulsestrikes a tilted planar slab at 45◦ immersed in vacuum: one in which the slab is made ofpositive-index glass and the other in which the slab is made of a hypothetical negative-index medium. As shown in Figure 4.9, global momentum and center-of-mass velocity areconserved using all five postulates regardless of whether the pulse undergoes a positive ornegative refraction in the slab. This example shows that the methods used in this thesis forforce density calculations are consistent with conservation laws in two dimensions. Directnumerical calculations offer an alternative to other approaches in multiple dimensions thatseparate force density terms in order to explicitly describe the effects of scattering and fieldgradients such as in [73].As shown in Figure 4.9 and 4.10, the five postulates make different predictions of theslab recoil and therefore, can potentially be distinguished by experiments. For the case ofa glass slab, the Minkowski postulates predict a momentary backward slab recoil, whereasall other postulates predict a forward slab recoil. This difference could be observed usingshort pulses in which recoil imbalances upon pulse entry and exit are separated in time.Measurements of pulse-induced recoil have yet to be performed, but according to thesesimulations could be a means to empirically validate a unique electrodynamic theory.40AB, EL, Chu,& AMP p0-p0MNAB, EL, & Chup0-p0MNAMP-p0 p0(a)(b)Pulse Slab Enclosure TotalFigure 4.9: Two-dimensional field-based implementation of the Balazs thought experiment.Non-degenerate momentum predictions made by the Minkowski, Abraham, Einstein-Laub,Chu, and Amperian postulates for closed systems in which finite pulses are obliquely inci-dent onto tilted impedance-matched slabs having (a) n =√3 (blue) and (b) n= −1.3+i0.07(yellow). The vector diagrams show instantaneous momentum of the pulse (red arrow), slab(blue arrow), enclosure (black arrow), and total system (green arrow) as the pulse propa-gates through the slabs. The total system momentum is conserved for all cases, as indicatedby the stationary green dots. Slab thicknesses (d = 5250 nm for the positive-index slaband d = 3500 nm for the negative-index slab) have been chosen so that the incident pulse(λo = 500 nm, τ = 2 fs) is fully immersed in the slab prior to exit. Dispersion of thenegative-index slab is described by a Drude model for the permittivity and permeabilitywith the following parameters: scattering rate Γ = 2250 THz, plasma frequency ωp = 1200THz, static permittivity ε∞ = 1.4, and static permeability µ∞ = 1.4.411 µmAB, EL, Chu, & AMP: SlabMN: SlabAB, EL, Chu, & AMP: PulseMN: PulseTotalAB, EL, & Chu: SlabMN: SlabAMP: SlabAB, EL, & Chu: PulseMN: PulseAMP: PulseTotal(a) (b)Figure 4.10: Non-degenerate center-of-mass displacement predictions made by theMinkowski (MN), Abraham (AB), Einstein-Laub (EL), Chu, and Amperian (AMP) postu-lates for closed systems in which finite pulses are obliquely incident onto tilted impedance-matched slabs having (a) n =√3 (blue) and (b) n˜ = −1.3 + i0.07 (yellow). The configura-tions are identical to those in Figure 4.9. Note that the total center-of-mass displacementremains fixed, as indicated by the stationary green dots.424.5 Radiation Pressure on Submerged MirrorsExperiments to study radiation pressure in the early 1900s by Lebedev [3], Nichols andHull [4, 5], and later by Bell and Green [74], proved that light impinging from a vacuumcould transfer momentum to macroscopic objects. It was not until the mid-1900s thatattempts were made to probe the momentum conferred by light impinging from dielec-tric media. In 1954, Jones and Richards [13] immersed a rhodium-coated silver mirror invarious dispersive liquids ranging in refractive index from 1.33 to 1.63, and observed itsdisplacement caused by tungsten lamp illumination. They concluded that light of a fixedintensity exerted a normalized pressure that is proportional to the refractive index of theliquids within an error margin of ±1.2%. This experiment was later refined in 1978 byJones and Leslie [14]. The lamp was replaced with a HeNe laser to increase light coherence,and the rhodium-coated silver mirror was replaced with a multi-layered dielectric mirror (ofunknown composition) to mitigate heating and thermal expansion. Improved experimentalprecision yielded confirmation, within an error margin of ±0.05%, of a linear dependence ofnormalized radiation pressure on the refractive index of the immersing fluids. Although theresults of the Jones-Richards and Jones-Leslie experiments are directly explained by howlight momentum of the Minkowski form scales with refractive index [25, 26, 27, 28], therehas been growing evidence that the experiments could be modelled using different electro-dynamic formulations [40, 42, 75, 76, 77, 78]. As early as 1978, Jones [40] proposed anAbraham model to explain momentum transfer to submerged mirrors based on mechanicalperturbations in the fluid which add to the electromagnetic contribution of momentum togive a Minkowski-like total momentum. Field-based descriptions of the experiments haveaccurately predicted the dependence of radiation pressure on refractive index using theEinstein-Laub postulates [77], the Minkowski and Chu postulates [28, 78], and the Lorentzforce density [75, 76]. These efforts, in turn, have been used as evidence to support theAbraham momentum density [77], the Minkowski momentum density [28, 78], or the arith-metic mean of the two [75, 76].Different routes can be taken to arrive at identical predictions of radiation pressure.This is illustrated by considering pulse interaction with a reflective mirror assumed to bea perfect electric conductor, a common approximation of the mirrors used in the Jones-Richards and Jones-Leslie experiments [77, 42, 28]. The pulse, mirror, and surroundingdielectric medium are within an enclosure that is initially stationary such that 1) dynamicmomentum must always achieve a zero sum and 2) the origins of the mirror momentumcan be visualized in time. We assume that the surrounding medium and enclosure areconnected. A pulse is created from the left end of the enclosure and bounces off the mirror43at normal incidence. As shown in Figure 4.11, the postulates predict dynamics that are allconsistent with global momentum conservation and indistinguishable except for one. TheMinkowski postulates, the lone outlier, describe the transfer of momentum upon reflectionas a two-body process involving just the pulse and mirror. The pulse is created withmomentum np0 and the system recoils with equal and opposite momentum −np0, wherep0 is the momentum of the pulse in vacuum. With the system momentum fixed at −np0over the remaining duration of the simulation, the pulse confers momentum 2np0 directlyto the mirror upon reflection. The other postulates (Abraham, Einstein-Laub, Amperian,and Chu) describe the transfer of momentum as a three-body process involving the pulse,surrounding fluid, and the mirror. The pulse is created with momentum p0/n and thesystem recoils with momentum −p0/n. When the pulse is reflected from the mirror, twomechanisms transfer momentum to the mirror. The mirror receives a parcel of momentum2p0/n due to reversal of the pulse direction and another parcel of momentum 2(n2−1)p0/npresent in the surrounding fluid, resulting in a final mirror momentum 2np0. Because thefinal mirror momentum for all postulates is the same, and differences in the mediatingmechanisms are likely immeasurable, the five postulates in this case are degenerate.44−np0np02np0PulseEnclosureMirrorTotal−p0/np0/n2np02p0(n2 − 1)/nt=1 fs t=16 fs t=28 fs t=37 fs0 50time (fs)(a)(b)MNAB, EL,Chu, & AMPnormalizedmomentumFigure 4.11: Momentum transfer between a pulse (λ0 = 632 nm, τ = 2 fs) in an enclosurecontaining a fluid (n = 1.6) and a submerged mirror approximated as a perfect electricconductor. The top panel shows simulation snap-shots of the electric field as the pulseis reflected from the mirror. The momentum of the pulse (red solid line), enclosure andsurrounding fluid (black dotted line), mirror (blue solid line), and total system (greensolid line) predicted by (a) the Minkowski postulates and (b) the Abraham, Einstein-Laub,Chu, and Amperian postulates. All five postulates predict the same final slab momentum of2np0. Note that the enclosure and surrounding fluid are connected such that the momentumassociated with both are lumped together.45Lorentz force calculations of radiation pressure on a submerged reflector imparting arbi-trary phase show that the acquired momentum can vary continuously with phase betweenan upper bound of the Minkowski momentum and a lower bound of the Abraham momen-tum [75, 76]. Pulse interaction with a perfect magnetic conductor – the compliment of aperfect electric conductor – illustrates how a half-cycle phase shift in the reflected electricfield can drastically alter radiation pressure predictions. Repeating the previous exampleusing a perfect magnetic conductor in Figure 4.12 shows that all postulates remain faithfulto global momentum conservation, but the degeneracy of the final mirror momentum is nowbroken. The Minkowski and Abraham postulates predict a final mirror momentum 2np0resulting from either two-body or three-body momentum transfer processes upon pulse re-flection. The Einstein-Laub, Chu, and Amperian postulates, on the other hand, predicta smaller final slab momentum of 2p0/n, which is acquired directly by the reversal of thepulse momentum without the accompanying system recoil observed for the case of a perfectelectric conductor.46−np0np02np0PulseEnclosureMirrorTotal−p0/np0/n2np02p0(n2 − 1)/n−p0/np0/n2np0t=1 fs t=16 fs t=28 fs t=37 fs0 50time (fs)(a)(b)(c)MNABEL, Chu,& AMPnormalizedmomentumFigure 4.12: Momentum transfer between a pulse (λ0 = 632 nm, τ = 2 fs)in an enclosurecontaining a fluid (n = 1.6) and a submerged mirror approximated as a perfect magneticconductor. The top panel shows simulation snap-shots of the electric field as the pulse isreflected from the mirror. The momentum of the pulse (red solid line), enclosure and thesurrounding fluid (black dotted line), mirror (blue solid line), and total system (green solidline) predicted by (a) the Minkowski postulates, (b) the Abraham postulates, and (c) theEinstein-Laub, Chu, and Amperian postulates. The Minkowski and Abraham postulatespredict a final slab momentum of 2np0, whereas the Einstein-Laub, Chu, and Amperianpostulates predict a final slab momentum of 2p0/n. Note that the enclosure and surroundingfluid are connected such that the momentum associated with both are lumped together.47Given that radiation pressure on submerged mirrors is not completely degenerate andcan be phase dependent, accurate modelling of the Jones-Richards and Jones-Leslie ex-periments must fully account for the specific configurations in which radiation pressure isgenerated. We use the five electrodynamic postulates in simulations that describe the twoexperiments to a degree of realism unmatched to date. The simulations incorporate realisticexperimental conditions such as broadband illumination, the skin effect in a metallic mirror,multiple reflections in a dielectric mirror, and frequency-dependent dispersion. To mimicpower normalization used in both experiments, the incident time-averaged power flux foreach case of immersing fluid is kept constant, a consideration whose importance has beenrecently pointed out [77].We begin by considering the Jones-Richards experiment in which incoherent, broadbandlight from a tungsten lamp illuminates a rhodium-coated silver mirror. Virtual replicationof the original experiment is limited by missing information, such as the spectral distribu-tion of the tungsten lamp, the complex permittivity of the rhodium-coated silver, and thethicknesses of rhodium and silver. In the absence of these details, we make the reasonableassumption that the mirror is made from either bulk rhodium or bulk silver, modelling eachmetal by fitting the Drude model to experimentally measured permittivity values found in[79] for rhodium and [80] for silver. The tungsten light source is modelled as a summationof incoherent plane waves with a Gaussian spectral distribution centered at the wavelength632 nm. Key assumptions of the analysis, including how the Drude model is fitted topermittivity values of rhodium and silver, and the assumed spectral distribution of thetungsten light source, are shown in Figure 4.13.48300 600 1200 1600 2000 2400−100−80−60−40−20020406080100λo (nm)ε′ r,ε′′ r400 100000.51Power(a.u)ε′Ag[80] ε′′Ag[80]ε′Ag ε′′Agε′Rh[79] ε′′Rh[79]ε′Rh ε′′RhFigure 4.13: Key assumptions used to model the Jones-Richards experiment. The real partof complex permittivity ε′ (solid lines) and the imaginary part of the complex permittivityε′′ (dotted lines) produced by the Drude model when fitted to tabulated permittivity valuesof rhodium (markers: triangle, cross) [79] or silver (markers: circle, square) [80]. The insetdepicts the power spectrum of a broadband incoherent light source assumed to emulatelight from the tungsten lamp used in the Jones-Richards experiment.The radiation pressure on a submerged metallic mirror predicted by different postulatescan change dramatically depending on the assumed properties of the mirror. As shown inFigure 4.14 (a), if the mirror is assumed to be a perfect electric conductor, all five postulatesidentically predict pressure that linearly scales with the fluid index, a result that matcheswithin the experimental error. Not surprisingly, the predicted pressure becomes inverselyproportional to the fluid index for the Einstein-Laub, Chu, and Amperian postulates whenthe perfect electric mirror is replaced by a perfect magnetic mirror. However, when wemodel the metallic mirror as a dispersive metal, the coincidence between theory and exper-iment seen before by assuming a perfect electric conductor is weakened. Incorporating thedispersion of the metal introduces two realistic effects: losses in the metal and a slight, butnon-zero, phase shift imparted by reflection. As shown in Figure 4.14 (b), the pressurespredicted by all five postulates exhibit departures from the Jones-Richards experimentaldata that become significant for higher fluid index values. This departure is even morepronounced if the mirror is rhodium as opposed to silver, due to the larger optical losses ofrhodium, which fall further from the perfect electric conductor approximation. For identical49simulation parameters, the Einstein-Laub, Chu, and Amperian postulates are generally inworse agreement to the experimental data than the Abraham and Minkowski postulates.However, without detailed knowledge of the exact composition of the rhodium-coated silvermirror, it is not possible to conclude whether any of the postulates accurately predict theJones-Richards results. Further radiation pressure measurements using well-characterizedmetallic mirrors are thus a necessary step toward the isolation of an electrodynamic theory.We next consider the Jones-Leslie experiment in which coherent radiation from a laserilluminates a submerged dielectric layered mirror, which results in pressure scaling linearlywith fluid index. Again, virtual replication of the original experiment is limited because thecomposition of the dielectric mirror (which was purchased from a vendor) was not reported.Here, we assume a dielectric mirror made from alternating, quarter-wavelength-thick lay-ers of MgF2 (n = 1.37) and ZnS (n = 2.35), a recipe reported by Jones and Leslie [14]to make some of the other mirrors used in their experiment. Interestingly, the orderingof the layers can significantly alter the radiation pressure predicted by some of the pos-tulates, as shown in Figure 4.14 (c). The Abraham and Minkowski postulates predict apressure proportional to n regardless of the layer ordering, but the other three postulatespredict either n-proportional or 1/n-proportional radiation pressure depending on whetherthe higher-index layer (ZnS) or the lower-index layer (MgF2) is the exterior layer in contactwith the surrounding fluid. The sensitivity of the radiation pressure to the layer orderingarises from the phase variations in the reflection from a dielectric mirror. That is, thephase of reflection from a dielectric mirror when capped by the higher-index layer is likethat from a perfect electric mirror. When capped by the lower-index layer, it is like thatfrom a perfect magnetic mirror. Without knowing the composition of the dielectric mirrorused in the Jones-Leslie experiment, we cannot rule out any of the postulates on the basisof the experimental accord, since all postulates are capable of predicting the n-proportionalradiation pressure in agreement with the experimental data. In light of this analysis, anessential experiment to test the hypothesis of phase-dependent radiation pressure is to mea-sure the radiation pressure difference on two submerged dielectric mirrors whose layers areordered differently but otherwise identical.500.60.81.01.21.41.61.82.02.2GRR〈Pz〉/〈Pz(n=1)〉All forms (PEC)AB & MN (PMC)EL, Chu, & AMP (PMC)Jones & Richards (1954)0.60.81.01.21.41.61.82.02.2GRR〈Pz〉/〈Pz(n=1)〉AB & MN (Ag)AB & MN (Rh)EL, Chu, & AMP (Ag)EL, Chu, & AMP (Rh)Jones & Richards (1954)1.0 1.1 1.2 1.3 1.4 1.5 1.6 1.70.60.81.01.21.41.61.82.02.22.42.6refractive index〈Pz〉/〈Pz(n=1)〉AB & MNEL, Chu, & AMP (ZnS/MgF2)EL, Chu, & AMP (MgF2/ZnS)Jones & Leslie (1978)1.0 1.1 1.2 1.3 1.4 1.5 1.6 1.70.60.81.01.21.41.61.82.02.22.42.6refractive index〈Pz〉/〈Pz(n=1)〉AB & MN (dispersive)AB & MNEL, Chu, & AMP (dispersive)EL, Chu, & AMPJones & Leslie (1978)(a) (b)(d)(c)Figure 4.14: Normalized time-averaged radiation pressure on submerged mirrors predictedby the Minkowski, Abraham, Einstein-Laub, Chu, and Amperian postulates. The greendots and error bars in (a) and (b) are data from the 1954 Jones-Richards experiment [13].The green dots in (c) and (d) are data from the 1978 Jones-Leslie experiment [14]. Notethat the error in the 1978 Jones-Leslie experiment is smaller than the dot size. The mirrors,immersed in fluids with refractive index n that varies from 1.0 to 1.7 in steps of 0.1, areilluminated at normal incidence. Pressure predictions are made for the following cases:(a) A perfect mirror (electric or magnetic) submerged in fluids without dispersion; (b) amirror made of dispersive silver or rhodium (whose permittivity is modeled using the Drudemodel fitted to tabulated data [80, 79]); (c) a multi-layered dielectric mirror composed ofalternating λ0/4-thick layers of MgF2 (refractive index 1.38) and ZnS (refractive index 2.35)submerged in fluids without dispersion fluids for configurations in which the outermostdielectric layer is MgF2 (MgF2/ZnS) or ZnS (ZnS/MgF2); (d) a multi-layered (MgF2/ZnS)dielectric mirror where the amount of dispersion in the fluids is set to ng ≈ 1.03n, whereng is the group refractive index and n is the phase refractive index (corresponding to themaximum upper bound of the fluid dispersion as reported by Jones and Leslie [14]). Themirrors in (a) and (b) are illuminated at normal incidence by an incoherent white lightsource with spectral intensity shown in the inset of Figure 4.13, which approximates atungsten lamp. The mirrors in (c) and (d) are illuminated at normal incidence at thewavelength λ0 = 632.8 nm.51An important conclusion of the Jones-Leslie experiment is that radiation pressure isdependent on the phase refractive index as opposed to the group refractive index. Wesimulate the pressure on dielectric mirrors (with a low-index capping layer) immersed influids that are modelled as dispersive dielectrics with variations in group and phase re-fractive index matching those reported by Jones and Leslie [14]. As shown in Figure 4.14(d), adding the effect of dispersion to the surrounding fluid medium has a negligible effecton the pressure trends for all postulates, which remain dependent on phase refractive index.Interpretations of radiation pressure on submerged mirrors are muddled by the presenceof force densities residing in the adjacent fluid. There is no consensus on the importanceof these contributions to the total radiation pressure on the mirror, as they have beenneglected by some [77] and incorporated by others [78, 81]. To illustrate their relativeimportance for each of the five postulates, Figures 4.15 (a) and (b) plot the time-averagedforce density exerted by a normally incident continuous wave onto a multi-layered dielectricmirror, capped by a low-index layer and a high-index layer, respectively. The Minkowskiand Abraham postulates predict a series of pressure spikes localized at the interfaces ofthe dielectric mirror but no force density in the fluid region. The Einstein-Laub, Chu,and Amperian postulates, on the other hand, predict oscillatory force density patterns thatspan across the fluid and dielectric mirror regions. Consistent with previous observations,the force density in the dielectric mirror is larger when the mirror is capped by the higher-index layer. The force density in the fluid region can add or subtract to the total radiationpressure on the submerged mirror. As shown in Figure 4.16, the pressure on the mirroroscillates between n and 1/n by accounting for the force density every quarter cycle intothe fluid, matching similar predictions made in [42, 75]. This adds an additional layerof complexity to the interpretation of the Jones-Leslie experiments, as the observed indexproportionality of radiation pressure measured by Jones and Leslie can be recovered by theEinstein-Laub, Chu, and Amperian postulates through selective force density integrationinto the fluid.520.00.20.40.60.81.0MgF2ZnSnormalizedforcedensity1.0 1.2 1.4 1.6 1.8 2.0 2.2 2.4-0.050-0.0250.0000.0250.0500.0750.100z (µm)normalizedforcedensity(a)ABMNEL & ChuAMP0.00.20.40.60.81.0normalizedforcedensity1.0 1.2 1.4 1.6 1.8 2.0 2.2 2.4-0.050-0.0250.0000.0250.0500.0750.100z (µm)normalizedforcedensity(b)Figure 4.15: Time-averaged force density distributions predicted by the Minkowski, Abra-ham, Einstein-Laub, Chu, and Amperian postulates, shown for the case of a continuous-wave beam (λ0 = 632.8 nm) normally incident onto a dielectric mirror composed of alternat-ing λ0/4-thick layers of MgF2 and ZnS submerged in a fluid of refractive index n = 1.6. Weconsider the configuration for which the dielectric mirror is capped by (a) the lower-indexlayer and (b) the higher-index layer. Degenerate time-averaged force density distributionspredicted by (top) the Abraham (solid blue) and Minkowski (dashed red) postulates and(bottom) the Einstein-Laub and Chu postulates (solid black) and the Amperian postulates(dashed orange). The former is dominated by strong surface forces localized at transi-tions from high to low refractive index regions, whereas the latter is dominated by weakeroscillatory body forces extending across the fluid and mirror regions. All force densitydistributions have been normalized to the peak instantaneous Abraham force density.53−5λo/4n −3λo/4n −λo/4n 01/nndielectricmirrorfluid, nelectrical length φz〈Pz〉/〈Pz(n=1)〉MgF2/ZnSZnS/MgF2Figure 4.16: Normalized radiation pressure on a submerged dielectric mirror predictedby the Einstein-Laub, Chu, and Amperian postulates as a function of the force densityintegration distance into the fluid for cases in which the dielectric mirror is capped byeither the lower-index layer (blue solid) or the higher-index layer (red dashed). Here, nrepresents the refractive index of the surrounding fluid medium. The circles denotes thenominal radiation pressure when the force density is integrated up to the boundary of themirror. The addition of body forces in the fluid causes the normalized radiation pressureto oscillate between an upper bound of n and a lower bound of 1/n.544.6 Radiation Pressure at Fluid InterfacesPerhaps the simplest configuration to probe the momentum of light is to observe itspassage across the threshold between vacuum and a dielectric medium. Changes in themomentum of light upon entering the medium should be offset by an-axis recoil dependenton whether momentum takes on Minkowski or Abraham forms [82, 83]. For the simple caseof a beam of intensity I normally incident from air onto a dielectric medium of refractiveindex n, the Minkowski form of momentum predicts a negative pressure (pull) of [82]P = −2Icn− 1n+ 1, (4.9)whereas the Abraham form of momentum predicts a positive pressure (push) ofP =2Icn− 1n+ 1. (4.10)In 1973, Ashkin and Dziedzic [15] performed an experiment to measure the longitudinalrecoil exerted by pulses entering a deformable dielectric medium. In their experiment, anair-water interface was illuminated with high-power (3 kW), 60-ns-long, 530-nm-wavelengthlaser pulses focused down to a spot size of radius 2 µm. Based on the far-field profile of thetransmitted pulse, they inferred that the pulse induced a momentary, micron-scale, upwardbulge on the water surface, overcoming both surface tension (' 7 mJ/m2 at room temper-ature) and gravity. The direction of the bulge, consistent for illumination downward fromthe air or upward from water, was initially interpreted by Ashkin and Dziedzic as evidenceof the on-axis recoil associated with the Minkowski form of light momentum [15]. Sub-sequently, Gordon [84] proposed an alternative explanation based on compressive lateralradiation forces squeezing the water into the highest-intensity region of the beam, a so-called “toothpaste” effect that forced water to bulge up into the air above. This viewpointhas gained traction [25, 85] but has yet to be proven. Recently, laser-induced deformationsof an air-water interface were shown to be accurately modelled without explicitly invokingcompressive forces [83, 86]. The model described interface deformation using the Navier-Stokes equation driven by a Minkowski-like on-axis recoil given by equation. 4.9, which isevidence that suggests a Minkowski form of momentum.Radiation-induced fluid deformations of much greater effect have been studied by Casneret al. [16, 17, 18, 19, 20]. In place of an air-water interface was a carefully designed liquid-liquid interface, where the composition of the liquids was tuned to achieve a surface tension104 times smaller than that of an air-water interface. Because the higher-index fluid sat55atop the lower-index fluid, deformation of the interface towards the lower-index fluid wouldbe directed along, rather than opposed to, gravity. It was shown that illumination of thefluid-fluid interface with even a moderately-powered ('150 mW) continuous-wave 632-nm-wavelength laser beam could cause millimeter-scale downward bulges into the lower-indexliquid, irrespectively of the illumination direction [19, 20]. The much larger bulge magni-tudes enabled detailed studies of the interface morphology and asymmetric deformationsdependent on illumination direction [87]. The similarities in the light-induced deforma-tions observed for air-water and liquid-liquid interfaces – both deform from high-index tolow-index regions – suggest a common physical origin, but a consolidated model rootedin the basic electrodynamic postulates has yet to be proposed. Such efforts, however, areimportant for resolving the roles of lateral and longitudinal forces and elucidating the link,if any, between fluid deformations and light momentum in matter.Realistic descriptions of radiation-induced fluid deformation must abandon the assump-tion of material rigidity, a significant departure from ideality that requires careful re-examination of the nature of recoil. Figure 4.17 shows snap-shots of the instantaneous forcedensity distribution described by all five postulates, in the case of a plane wave impingingat normal incidence onto a dielectric from free space. The force density distributions arecharacterized by localized surface pressure (Minkowski), undulating body forces (Einstein-Laub, Chu, and Amperian), or a combination of both (Abraham). Under the assumptionthat the medium is rigid, momentary imbalances in the force density distribution extendingacross the entire medium can result in a net recoil. However, under the realistic assumptionthat the medium is deformable, such recoil mechanisms are impossible as it requires instan-taneous communication of force density across finite extents. Instead, each element of themedium is locally driven by a time-averaged force density, which, as shown in Figure 4.17,can be substantially different from the instantaneous force density. Time-averaging gener-ally diminishes the undulations in the instantaneous force density distributions, so much sofor the Abraham force density that it becomes indistinguishable from the Minkowski forcedensity. In this case, models of recoil based on the Minkowski and Abraham time-averagedforce densities predict identical backwards recoil driven entirely by surface pressure. Thiscontrasts to classical models of recoil based simply on changes in light momentum at thresh-old of a dielectric, which predict opposing recoils for light momentum of Minkowski andAbraham forms [82, 83].Both electrodynamics and fluid dynamics are needed to completely understand howfluids deform under radiation pressure. Since the time scales of fluid deformation aremuch greater than the times scales of electromagnetic oscillation (particularly at visible56fpk0f z(z,t)ABMNEL & ChuAMPfpk0z (µm)〈fz(z)〉z (µm)(a) (b)Figure 4.17: Instantaneous (top row) and time-averaged (bottom row) electromagnetic forcedensity distribution exerted by a continuous-wave beam (λ0 = 632.8 nm) normally incidentfrom air onto water (n = 1.33) predicted by (a) the Abraham (solid blue) and Minkowski(dashed red) postulates and (b) the Einstein-Laub and Chu postulates (solid black) andthe Amperian postulates (dashed orange).frequencies), the physics of radiation-induced fluid deformation can be captured withoutsolving the electrodynamic and fluid dynamic equations simultaneously. We use a cou-pled electrodynamic-hydrodynamic simulator in which Maxwell’s equations and the Navier-Stokes equations are solved on respective grids of fine and course spatio-temporal resolution.The electrodynamic simulator calculates the time-averaged force density, which is portedto the hydrodynamic simulator to calculate an updated fluid geometry. We examine twoconfigurations: pulsed excitation of an air-water interface like in the Ashkin-Dziedzic exper-iment [15] and continuous-wave visible excitation of a low-tension fluid-fluid interface likein the experiments by Casner et al. [19, 20]. To emulate impulsive excitation, we drive thehydrodynamic simulator at a single coarse time step with the time-averaged force densityprofile calculated for the initial fluid configuration, which is a valid procedure only becausethe coarse time step is several orders of magnitude longer than the pulse duration. Toemulate continuous-wave excitation, we cyclically compute the time-averaged force densityand resulting fluid geometry each time the interface has deformed by at least one coarsespace step. The simulations assume a simplified two-dimensional configuration which doesnot fully describe the tightly focused beams used in the experiments [15, 19, 20], but hasa sufficient degree of complexity to illustrate lateral forces associated with electrostrictionand field gradients.57Figure 4.18 illustrates the initial time-average force distribution applied to an impulsively-excited air-water interface and the ensuing interface deformation and fluid velocity. Threeunique sets of bulge dynamics emerge. For both transverse-magnetic (TM) and transverse-electric (TE) polarizations, only the Abraham, Minkowski, and Einstein-Laub postulatespredict an upward bulge for illumination from air or water, qualitatively consistent withthe results of the Ashkin-Dziedzic experiment [15]. The Amperian and Chu postulates, onthe other hand, predict upward bulges for TE polarization, which are smaller than thosepredicted by the other postulates, and downward bulges for TM polarization. The lattercase provides further proof of the empirical invalidity of the Amperian and Chu postulates.This conclusion is consistent with those made by Mansuripur et al. [88] after comparingthe effect imparted by a beam to a dielectric medium using the Einstein-Laub and Lorentzforce densities.58AB & MNELAMP & Chu: TMAMP & Chu: TEBeam directed down Beam directed up〈~fe〉 ~vFigure 4.18: Two-dimensional fluid dynamic simulations representative of the Ashkin-Dziedzic experiment [15]. The lower fluid (blue) represents water (n=1.33) and the upperfluid (white) represents air (n=1). The air-water interface is excited by a 60-ns-long, 530-nm-wavelength pulse at normal incidence. The time-averaged force density distributionsexerted by the pulse are calculated using the Minkowski, Abraham, Einstein-Laub, Chu,and Amperian postulates (magenta arrows), which are shown next to the resulting velocityfield (blue arrows) of the deformed interface. We set the surface tension of the water-airinterface to σ = 0.7 mN/m and the viscosity of water to η = 1 mPa·s. The Abraham,Minkowski, and Einstein-Laub postulates predict an upward bulge for TE and TM polar-izations for both illumination directions. The Amperian postulates predict an upward ordownward bulge for TM and TE polarizations, respectively. Mesh refinement studies havebeen performed to ensure convergence of the simulated fluid behaviour.59The Abraham, Minkowski, and Einstein-Laub postulates offer two physical mechanismsto explain the appearance of an upward bulge: it can arise from a longitudinal upwardforce density localized at the interface (Minkowski and Abraham) or a lateral compres-sive force density applied to the bulk of water (Einstein-Laub). The former is consistentwith the Minkowski surface recoil described by equation 4.9, and the latter describes thetoothpaste effect suggested by Gordon [84]. The bulge shape and size resulting from ei-ther mechanism are nearly identical and indistinguishable with respect to the precision ofthe Ashkin-Dziedzic experiment [15], although repetition of the experiment with fluids ofdifferent viscosities and mass densities may lead to more diverse bulge shapes that betterdistinguish surface and body effects. For the Abraham, Minkowski, and Einstein-Laubpostulates, illumination from water yields larger bulges than illumination from air, a gen-eral consequence of the larger field magnitudes and larger force densities in water whenit is used as the incidence region. The interplay between force density distributions andinterface deformation can be complex. For example, the Einstein-Laub postulates predictdownward body forces when illuminated from water, which are overwhelmed by an upwardsurface force and overall compressive forces that conspire to generate an upward bulge.The bulge dynamics observed by Casner et al. [19, 20] for a fluid-fluid interface havebeen distinguished by two regimes: a low-power regime in which the bulge sizes are smalland scale linearly with power and a high-power regime in which the bulge exhibits largeinterface instabilities. Recent investigations have suggested that the instabilities in thehigh-power regime are likely caused by surface heating [87], which cannot be treated byour simulations. We therefore restrict our analysis to the low-power regime. Figure 4.19shows two sets of simulated bulge and force distribution dynamics: one predicted by theAbraham or Minkowski postulates and the other by the Einstein-Laub postulates. Thebulge dynamics of the fluid-fluid interface are consistent with the dynamics of the air-waterinterface, but reversed in direction since the optically denser region now resides above theinterface. The downward bulge for both illumination directions and the slightly larger bulgesizes for downward illumination are both consistent with the observations [19, 20]. The bulgecan originate from a downward surface pressure (Minkowski and Abraham) or a pinchingeffect distributed throughout the bulk (Einstein-Laub). Bulges due to surface pressureexhibit larger asymmetries with respect to illumination direction, a feature potentiallyuseful for empirical validation of the two mechanisms of bulge formation.60AB & MN: ↓AB & MN: ↑EL: ↓EL: ↑Figure 4.19: Two-dimensional fluid dynamic simulations representative of the experimentsby Casner et al. [19, 20]. These experiments used a micellar phase of microemulsions madefrom a liquid mixture of mass composition 9% water, 4% sodium dodecyl sulfate, 70%toluene, and 17% n-butanol. Below the critical temperature of 35 C◦, the mixture separatesinto two phases of different micellar concentrations. The two phases are represented in thesimulations as two continuous fluid media in which the upper fluid (beige) has index n= 1.49 and mass density ρ = 853 kg/m3 and the lower fluid (blue) has index n = 1.46and mass density ρ = 1017 kg/m3 [89]. Near the critical temperature, the fluids havethe same viscosity of η = 1.3 mPa·s and an extremely low surface tension of 0.1 µN/m.The fluid mixture is excited by a continuous-wave TM-polarized beam (λ0 = 532 nm)with a Gaussian profile. The time-averaged electromagnetic force density distributions(magenta arrows) calculated using the Minkowski, Abraham, and Einstein-Laub postulatesfor downward (↓) and upward (↑) illumination are shown for three instances during theevolution of a radiation-pressure-driven surface bulge. The three postulates all predictdownward bulge formation from the high-index region into the low-index region, consistentwith the observations by Casner et al.[89]. Mesh refinement studies have been performedto ensure convergence of the simulated fluid behaviour.614.7 Optical TweezingThe most common application of radiation pressure is for moving small objects in low-friction environments, a practice known as optical tweezing. Popular examples are foundin biology, such as trapping of bacteria for selective destruction [90] and measurement ofthe forces involved in RNA transcription [91]. Optical tweezers are generally considered tooperate using two different kinds of forces: gradient forces, and scattering forces. Gradientforces refer to a particle’s tendency to move towards the location of highest beam intensity,and scattering forces are related to forces produced during reflection and transmission atthe interface between two media. Gradient forces can be thought of as primarily dependenton the spatial derivative of the electromagnetic energy density, where as scattering forcescan exist in the absence spatial gradients. Interestingly, one way to understand gradientforces is based on energy minimization. Recall in Section 2.3.1 where the force on a piece ofsteel was found by considering a force equal to the work needed to overcome the change inmagnetic energy inside the iron relative to that just outside it. Gradient forces on dielectricparticles can be understood in a similar way, where the energy of an optical tweezing systemis minimized when a dielectric particle is at the center of highest beam intensity [92]. Foroptical tweezers, the magnitude of this gradient force relative to the scattering forces alsohas to be considered, an interesting topic that is outside of the scope of this thesis. Anavenue of potential further work is to use the numerical methods in this thesis to quantifythe different magnitudes of scattering and gradient forces for various optical tweezing con-figurations. This would be an important step towards a complete theory of optical tweezers,which has yet to be developed [93, 48].Barton et al. [94] pioneered the approach of using Mie theory to express the fields scat-tered from a spherical particle in conjunction with the electromagnetic stress tensor tocalculate the force on that particle. The original expressions derived by this approach [94]used the Minkowski form of the stress tensor, whose use was based on the work by [25]and [27]. The decision to use the Minkowski stress tensor over the alternatives has littleconsequence. As shown previously, the Minkowski form is identical to the Abraham pre-dictions for all time-averaged cases considered so far. Pfeifer et al. [95] have also arguedthat the Minkowski and Abraham postulates describe optical tweezers equivalently, differ-ing only with respect to ease of implementation. This seems to contradict calculations byMansuripur et al. [88] of the force density within illuminated dielectric particles, whichsuggest that various postulates can be distinguished by differences in their distributions offorce density.62To explore how different postulates describe optical tweezer effects, Figure 4.20 (a)depicts simulations that use all five postulates to model the interaction between a focusedcontinuous-wave laser beam with a small submerged dielectric cylinder placed at the edgeof the beam’s focus. For both TM and TE polarizations, only the Minkowski, Abraham,and Einstein-Laub postulates are consistent with optical trapping: the beam exerts a time-averaged force that drags the particle into the focus and pushes it along the direction ofbeam propagation [96]. The Amperian and Chu postulates, on the other hand, predictparticle trapping for TE polarization and particle expulsion for TM polarization. Forcespredicted by the Einstein-Laub postulates also exhibit polarization sensitivity, but to amuch lesser extent. The reason for this difference is because the Einstein-Laub, Chu, andAmperian forms have dominant body force terms that are dependent on the interferencepattern within the particle. In contrast, the Abraham and Minkowski forms are dominatedby momentum exchange at the interface between two media, as demonstrated previouslyin Figure 4.15.63All Forms: TEAB & MN: TMAMP & Chu: TMEL: TM1 µmpolysterene beadFigure 4.20: Simulated radiation pressure on a polystyrene bead (green) placed at the edgeof the focal point of a continuous-wave beam (λ0 = 532 nm) in a background medium ofwater (n = 1.33). The polystyrene bead (n = 1.58) has a diameter of 820 nm and the forceacting on it has been determined under the rigid body assumption by integrating the forcedensity distribution within the bead calculated using the Minkowski, Abraham, Einstein-Laub, Chu, and Amperian postulates. The bead is shown immersed in a backgroundelectromagnetic field calculated using the FDTD method.Experiments by Ashkin [96] have shown that air bubbles near a focused beam are re-pelled from the focus and pushed along the direction of beam propagation. Whereas theactuation of a dielectric particle is driven by optical forces acting on or within the particle,the actuation of an air bubble is driven by optical forces acting on the fluid region sur-rounding the bubble. To compare experimental observations with theory, we simulate thedynamics of an air bubble near a focused beam using the electrodynamic-hydrodynamicsimulator. As shown in Figure 4.21 (b), the five postulates all predict that the air bubbleis displaced away from the beam focus and along the direction of beam propagation, whichis qualitatively consistent with past observations by Ashkin [96]. Predictions of bubbleexpulsion differ in terms of rate, shape, and trajectory. Similar to previous examples, thefive postulates applied to the case of TM polarization yield three unique sets of dynamics:one degenerate for Abraham and Minkowski postulates, another for the Einstein-Laub pos-tulates, and a third degenerate for Amperian and Chu postulates. With all things equal,lateral bubble expulsion is predicted to be fastest for the Einstein-Laub postulates andslowest for the Amperian and Chu postulates, discrepancies which can potentially be vali-dated by measurements.64In recent years, there has been increased interest in using light beams to “pull” on aparticle such that the projection of its displacement onto the incident beam direction hasa negative sign. This has been referred to as “tractor beaming” [97, 98, 99, 100, 101, 102,103, 104]. For example, Kajorndejnukul et al. [103] showed that oil droplets confined at theinterface of air and water can be pulled along the interface and into the direction of a lightlyfocused beam obliquely incident onto the droplets from the air region. This observationwas interpreted to arise from an increase in the forward light momentum inside the droplet,propelling the droplet in the backward direction. Since the refractive index of the dropletis higher than that of air or water, an increase of light momentum in the droplet was at-tributed to a Minkowski form of momentum. Simulations of these experiments suggest thatalternative explanations are feasible. As shown in Figure 4.22, the Minkowski, Abraham,and Einstein-Laub postulates all predict that an oil droplet resting at the interface of airand water can be pulled into the direction of an obliquely incident laser beam, despite dif-ferences in how they predict light momentum to scale with refractive index. The Amperianand Chu postulates, on the other hand, predict a push.xyAB & MNxyELxyAMP & ChuFigure 4.21: Simulated dynamics of an air bubble placed at the edge of the focal point of aTM-polarized, continuous-wave beam (λ0 = 532 nm) in a background medium of water (n= 1.33, σ = 0.7 mN/m, η = 1 mPa·s). The air bubble (n = 1.00) has a diameter of 820 nm.The evolution of the deformable bubble in a background electromagnetic field is modelled intwo-dimensions by solving both electrodynamic and fluid dynamic equations. The simulatedbubble dynamics predicted by (top row) the Minkowski and Abraham postulates, (centerrow) the Einstein-Laub postulates, and (bottom row) the Chu and Amperian postulates.The dashed horizontal line indicates the axis of the focused beam. In all cases, the bubbleis pushed away from the beam axis and along the direction of beam propagation, albeitwith different bubble shapes and trajectories. Overall, the simulated behaviour using allpostulates is qualitatively consistent with observations by Ashkin [96].65θiAB & MNELAMP & Chu1 µmFigure 4.22: Tractor beaming of an oil droplet (n = 1.42) resting at the interface betweenair (n = 1.0) and water (n = 1.33). The length, height, and radius of curvature of thedroplet are 6.4 µm, 2 µm, and 6 µm, respectively, which are comparable to the dimensionsof the droplets studied in [99]. A TM-polarized plane-wave (λ0 = 532 nm) is incident ontothe droplet at an angle of incidence θi = 50◦. Resulting forces predicted by the Minkowskiand Abraham postulates (blue arrow), the Einstein-Laub postulates (dashed black arrow),and the Amperian and Chu postulates (magenta arrow). The reflected fields bouncing offthe droplet and propagating behind the diagonal source line have been removed for clarity.In another example of tractor beaming, Brzobohaty´ et al. [104] showed that two obliquely-incident plane waves can exert pulling forces on a dielectric particle placed in the regionof beam interference. This observation was modeled by Barton et al. using the Minkowskistress tensor [94]. Figure 4.23 (a) and (b) shows simulations of a dielectric particle illu-minated by two obliquely incident plane waves that work together to pull on the particle.Due to the highly symmetric illumination conditions that establish standing waves, all fivepostulates make identical force predictions for both polarizations. As shown in Figure 4.23(c), the angle- and polarization-dependence of the normalized lateral force match the resultsfrom calculations presented in [104], all of which indicate a tractor-beaming effect for TEpolarization and at steep angles of incidence. Thus, recent studies of tractor beaming do notoffer experimental proof to pin down a single electrodynamic theory. Also demonstrated inthis section is that a large portion of optical tweezer phenomena can be calculated using twodimensional finite difference techniques as opposed to complicated Mie theory approachescommonly implemented [104].6605.586000e+00θi~k1~k21 µmxy05.586000e+001 µmFpull68 70 72 74 76 78 80 82 84 86 88 90 920.00.51.01.52.0θi◦normalizedforce(c)(a) (b)All forms: TE-FDTDAll forms: TM-FDTDRef. [104]: TERef. [104]: TMFigure 4.23: Tractor beaming of a dielectric particle using two obliquely incident planewaves. Simulation snapshots of two TE-polarized plane waves (λ0 = 532 nm) incident atθi = 70◦ onto a polystyrene cylinder (d = 820 nm, n = 1.58) immersed in water (n = 1.33),shown (a) before the plane waves strike the bead and (b) at steady state when the planewaves collude to exert a pulling force in the −x direction. (c) Force on the bead calculatedby FDTD using the Minkowski, Abraham, Einstein-Laub, Chu, and Amperian postulatesfor both TE (blue circles) and TM (orange circles) polarizations as a function of incidentangle θi. The forces have been normalized to the force magnitude at θi = 70◦. Under thehighly symmetric illumination conditions studied here, all five postulates make degenerateforce predictions. Pulling forces are only achieved for TE polarization within the range76◦ < θi < 90◦. Excellent agreement is observed with the three-dimensional Mie-theorycalculations (solid lines) by Brzobohaty´ et al. [104], despite the reduced dimensionality ofthe simulations.674.8 Radiation Pressure in Negative Index MaterialsIn 1968 V. Veselago proposed that materials with negative values of µ and could beused to produce optical effects not possible with naturally occurring materials [105]. Mostof his predicted effects, such as negative refraction and a reversed doppler effect, have beenrealized [101, 106]. One that remains to be tested is his prediction that a mirror immersedin a negative index medium would be pulled towards the light source. Veselago based hisprediction of negative radiation pressure on a description of momentum density derived byRytov [107] (valid for a dispersive medium under the slowly varying envelop approximation)given by~GR =εµc2~E × ~H +~k2[∂ε∂ω~E2 +∂µ∂ω~H2], (4.11)where ~k is the wave vector associated with the carrier phase velocity of the wave. If electro-magnetic momentum density is proportional to phase velocity, then plane waves in negative-index media would carry momentum opposite to the direction of energy propagation. As acorollary, it was proposed that reversal of wave momentum would cause illuminated objectsimmersed in a left-handed medium to be pulled towards the light source. It will be shownhere that the Abraham and Minkowski postulates do indeed predict a pull on a mirror forlight incident from a negative index medium. However, this pull is likely not observable asit is smaller than the push felt by the surrounding medium during the interaction.It has not been established that the momentum of an electromagnetic wave in a left-handed medium reverses direction relative to energy propagation, as the validity of equa-tion 4.11 has been challenged [108]. We will use simulation to compare equation 4.11 withthe five electrodynamic postulate within Table 1.1 and see if we may consider electromag-netic momentum as negative within left handed media. To calculate the magnitude of thecarrier wave vector ~k, we use the spatial fourier transform of the electric field profile anduse the wave vector of peak intensity. Figure 4.24 tracks the electromagnetic momentumof a pulse entering a negative index medium from free space. Shown in Figure 4.24 (a) to(d), all postulates remain positive except for the Amperian form of momentum. This isbecause the Amperian definition of momentum is ~E × ~B, and a property of ideal negativeindex materials is that ~E and ~B are 180◦ out of phase.68-po0poAbraham-po0poMinkowski-po0poEinstein-Laub, & Chu-po0poAmperian5 10 15 20 25 30 35 40 45 50 55-po0poTime (fs)VeselagoPulse SlabEnclosure Total0.75 fs 14 fs 33 fs 45 fsFigure 4.24: A perfectly absorbing enclosure is partially filled with a negative-index slab.The enclosure produces a pulse and the electromagnetic momentum density of the pulseis tracked as it is born from within the left end of the enclosure and then enters thenegative index slab. Momentum predictions shown are the: (a) Abraham, (b) Minkowski,(c) Einstein-Laub & Chu, (d) Amperian, and (e) Veselago forms. All forms predict adecrease of the pulse electromagnetic momentum, but only the Veselago and Amperianforms are negative. The refractive index is set to n = −1.2 and impedance matched usingthe Drude model. The Drude parameters are ε∞ = µ∞ =√1.2, ωp=4.53× 1015, andγe ≈ 0, making the medium approximately lossless.69Now we explore the prediction of negative radiation pressure, Figure 4.25 shows a pulsecreated within the same negative index medium as in Figure 4.24, but is now reflected froma mirror, replicating the configuration of negative radiation pressure proposed by Vese-lago. For both the Abraham and Minkowski cases, indeed the mirror is pulled as Veselagopredicted. Unfortunately, this interesting effect is likely unobservable as the mechanicalmomentum imparted to the surrounding negative index medium after reflection is largerthan the momentum imparted to the mirror. Each form describes reflection as a three bodyprocess involving the mirror, surrounding medium, and pulse momentum. If we considerthe net momentum imparted to both the mirror and surrounding medium, the net effectwould be positive. Again we see that the force on the mirror cannot be inferred from theperspective of assuming the momentum from Figure 4.24 reflects and the mirror receivesmomentum equal to this change. This is similar to what was shown in the Section 4.5,where the mechanical component of momentum within the surrounding medium for eachform other than Minkowski was needed to calculate the total force on the mirror.70−2po02po 3.2po−2.6poAbraham−2po02po3.8po−2.6poMinkowski−po0po−0.4po1poEinstein-Laub & Chu10 15 20 25 30 35 40 45 50 55 60 65 70−5po05po−6.9po7.1poTime (fs)Amperian6 fs 25 fs 45 fs 60 fsPulse MirrorFluid (n < 0) TotalFigure 4.25: A perfectly absorbing enclosure is filled with negative index fluid and anAg mirror. The enclosure produces a pulse and the momentum of the pulse, fluid, andmirror is tracked. Momentum predictions shown are the: (a) Abraham, (b) Minkowski, (c)Einstein-Laub & Chu, and (d) Amperian forms. Interestingly, the Abraham and Minkowskipostulates do indeed predict a pull on the mirror, but this pull is smaller in magnitude thanthe net push received by the fluid after reflection. The Einstein-Laub and Amperian formsdiffer in that the mirror is pushed and the fluid is pulled. The negative index materialparameters are the same in as in figure 4.24.71We create an analog to the Jones and Leslie section by simulating the pressure depen-dence of a mirror submerged in media with varying refractive index spanning from positiveto negative values. The results are shown below in Figure 4.26. For the positive index caseof a PEC, all forms predicted the same trend of an increase in pressure with the refrac-tive index. However, for the negative index regime, we see very different behavior. TheAbraham and Minkowski postulates predict an identical increase in negative pressure pro-portional to the index. The Einstein-Laub and Chu forms interestingly no longer agree withany other form, and produce essentially no change in pressure within the negative indexregime. The Amperian form predicts a push that grows non-linearly as the surroundingmedium becomes more negative.−1.6 −1.4 −1.2 −1−1.8−1.1−0.40.311.72.43.13.84.5〈Pz〉/〈Pz(n=1)〉1 1.2 1.4 1.6−1.8−1.1−0.40.311.72.43.13.84.5refractive index n/ // /All Forms (+) Jones&Richards (1954)AMP (−) EL & Chu (−)AB & MN (−)Figure 4.26: Demonstration of the momentum a mirror experiences with respect to itssurrounding refractive index. All forms predict degenerate results for a refractive indexn > 1, but differences arise for n < 1. The Abraham and Minkowski postulates predictdegenerate trends of an increasing pull force as the index of the surrounding fluid becomesmore negative. The Einstein-Laub and Chu forms predict no appreciable change in forceas the index of the surrounding fluid becomes more negative. The Amperian form predictsan exponential increase in pressure as the refractive index becomes more negative.724.8.1 Optical Tweezing With a Negative Index BeadTo further demonstrate how radiation pressure effects are different when negative indexmedia are involved, we can look at how optical tweezing effects change when the objectbeing “tweezed” has a negative refractive index. As shown in Section 4.7, a particle will drifttowards the location of highest field intensity. Figure 4.27 (a) shows simulation results fora Gaussian beam illuminating a particle located ≈ 1µm from the beam axis and submergedin a medium with refractive index n = 1.33. Figure 4.27 (b) shows another simulationwith identical parameters as (a), except the bead has a negative index of refraction. Thenegative index bead force results are different for all forms except the Amperian form, whichstill predicts the same push away from the beam center as in the positive index case. Theother forms of Abraham, Minkowski, and Einstein-Laub, which previously predicted a pulltowards the beam center, now predict a push. The Einstein-Laub form still predicts a pulltowards the beam center however its magnitude relative to the other forms is quite reduced.AB & MNAMP & ChuELn > 01 µm(a)AB & MNAMPChuELn < 01 µm(b)ChuELFigure 4.27: Optical tweezing of an impedance matched negative index bead under TMillumination. The bead is placed at the edge of the focal point of a continuous-wave beam(λ0 = 532 nm) in a background medium of water (n = 1.33). The bead has a diameter of820 nm and the net force is solved for the bead for two cases, positive index in (a) with(n = 1.58) and negative index in (b) with (n = −1.58) the force acting on it has beendetermined under the rigid body assumption by integrating the force density distributionwithin the bead calculated using the Minkowski (MN), Abraham (AB), Einstein-Laub (EL),Chu, and Amperian (AMP) postulates. It can be seen that the positive index results in(a) from [1] are significantly different for the negative index bead in (b). The Abrahamand Minkowski, Amperian, and Chu postulates predict that a negative index bead in wateris pushed away from the beam center, while the Einstein-Laub form predicts a pull intothe beam focus. The negative index bead’s drude parameters are set to ε∞ = µ∞ = 1.1,ωe = ωm = 6.21× 1015, and γe = γm ≈ 0.This hypothetical optical tweezing experiment is not experimentally feasible with presenttechnology. However, as technology advances to possibly include situations such as Fig-ure 4.27, this thought experiment reaffirms our current understanding of electrodynamic73modeling needs to also advance. This was another example of how the five forms listedin Table 1.1 are experimentally distinguishable and enables further work to be done inunderstanding the behavior of each form in the negative index regime.74Chapter 5: ConclusionThis thesis presents a consolidated analysis of radiation pressure observations usingclassical electrodynamics described by five different formulations. The analysis includedtheoretical considerations of center-of-mass translation presented by Balazs [29], radiationpressure measurements on mirrors [13, 14], observations of radiation-induced deformationof fluid interfaces [15, 16, 17, 18, 19, 20], optical tweezers, [96, 103, 104], and thoughtexperiments involving radiation pressure in negative index materials [105]. Out of thefive electrodynamic theories considered, namely the Abraham, Minkowski, Einstein-Laub,Amperian, and Chu forms, at least two were shown to reproduce each observable considered.Each electromagnetic formulation was tested against their respective ability to conserveenergy, momentum, and center-of-mass velocity. Within this thesis, the results obtainedwere:• The symmetric form of the Minkowski and Amperian formulations can conservecenter-of-mass velocity, but both require a surface energy term whose existence mustbe rationalized. A simple experiment was carried out where this surface energy waspredicted to present itself as a measurable temperature difference across the entranceand exit faces of an illuminated waveguide. No temperature difference was observed.• Section 4.5 simulated past experiments done to measure the force on metallic and di-electric mirrors while submerged in liquids with different values of refractive index [13,14]. For the 1978 dielectric mirror experiments, it was shown that the Einstein-Laub,Chu, and Amperian forms can predict a decrease in pressure with the surroundingrefractive index depending on the order of the dielectric layers, the Abraham andMinkowski forms maintain the same behavior for both possible orderings. The layerordering was never recorded for the original experiment, this means a second experi-ment could be done using both possible layer orderings.• The two dimensional analysis in Section 4.6 demonstrated that the Chu and Amperianforms could be differentiated by repeating the Ashkin or Casner experiments andcomparing results for TM and TE illumination.75• The results of two optical tweezing experiments were reproduced, the tendency ofa particle being tweezed to move towards the beam’s focus [96], and the ability oftwo interfering beams to pull on a particle [104]. Results agreed with [95] that boththe Abraham and Minkowski forms are sufficient to predict the forces involved inoptical tweezing. Polarization dependence was found for the Einstein-Laub, Chu, andAmperian forms that may be testable.• Lastly, in Section 4.8 we considered a thought experiment done by V. Veselago whopredicted electromagnetic momentum in negative index media could be considerednegative, and that a mirror immersed in a negative index medium would be pulledtowards the light source. It was shown that for an ideal case only the Amperianform predicts an electromagnetic momentum density that is negative. For his secondprediction we demonstrated that the Abraham and Minkowski forms do indeed pullon a mirror, but the pull force is less than the push that exists outside of the mirror inthe surrounding medium. As well, when a negative-index analogue of the Jones andLeslie experiments was re-created in Section 4.8 very distinct differences were seenbetween the five forms. It was also shown that optical tweezing effects reversed froma pull towards the beam focus to a push for the Abraham and Minkowski forms if theobject being moved was negative index.Rigorous simulations of five electrodynamic models demonstrated several situationswhere each model yielded different predictions. This large catalog of results agree withthe prevalent argument that the definition of electrodynamic theory is not arbitrary. Dif-ferences in both theoretical and experimental predictions of each form were shown andplausible experiments have been suggested to isolate a single true electrodynamic theory.76Appendix A: Speeding up Matlab CodeFigure HandlesMatlab is used extensively in research applications. When used for simulation, real timeplotting is common to help the user visualize their results. It is not immediately obviousto the typical Matlab user that their code can run much faster depending on how datais plotted in a figure window. Updating the handle on a plot can speed up your codeexponentially. An example of a slow plotting method and a fast plotting method are shownbelow in Figure 5.1, speedup will vary between computers but was ≈10 times faster for thecomputer used in this example. If your code uses many surface plots in two dimensions,it would not be unusual to see speedups of 500% by using handles instead of recalling the”figure” function every time a plot is to be updated.771 %% Slow Plotting2 f=60; % Frequency [ Hz ]3 T=1/f ; % Period [ s ]4 dt=T /200 ; % Time step [ s ]56 t=[0: dt : 4∗ T ] ; % Time axis [ s ]7 Nt=length ( t ) ; % # of points [#]8 E=zeros (1 , Nt ) ; % Electric field [ V/m ]910 f o r n=1:Nt1112 E ( n )=sin (2∗ pi∗f∗t ( n ) ) ;13 figure (1 )14 plot ( t ( 1 : Nt ) , E ( 1 : Nt ) ) ;15 xlabel ( ‘ time ‘ )16 ylabel ( ’E− f i e l d ’ )1718 end1 %% Fast Plotting2 f=60; % Frequency [ Hz ]3 T=1/f ; % Period [ s ]4 dt=T /200 ; % Time step [ s ]56 t=[0: dt : 4∗ T ] ; % Time axis [ s ]7 Nt=length ( t ) ; % # of points [#]8 E=zeros (1 , Nt ) ; % Electric field [ V/m ]910 figure (1 )11 plot_handle=plot (t , E ) ;12 xlabel ( ’ time ’ )13 ylabel ( ’E− f i e l d ’ )1415 f o r n=1:Nt16 E ( n )=sin (2∗ pi∗f∗t ( n ) ) ;17 set ( plot_handle , ’XDATA’ ,t , ’YDATA’ , E ) ;18 endFigure 5.1: A minimum working example (MWE) to demonstrate how to use handles toupdate a figure as opposed to re-calling the “figure” function. The code above will producea figure of a simple sine wave.78VectorizationWhen possible, a second method is to replace any loops in the code with a vector ,Matlab is more efficient at dealing with matricies as vectors rather than loops. Figure 5.2below shows Matlab code that will plot a 2D surface of a plane wave varying in time, thecode takes ≈ 10 seconds to run using for loops, and ≈ 4 seconds when vectorized as shown.1 f=60; % Frequency [ Hz ]2 T=1/f ; % Period [ s ]3 lambda=10; % Wavelength [ m ]4 k=2∗pi/lambda ; % Wavevector [ rad/m ]5 dt=T /90 ; % Time step [ s ]6 dx=lambda /30 ; % Space step [ m ]7 dy=lambda /30 ; % Space step [ m ]89 t=[0: dt : 4∗ T ] ; % t axis [ s ]10 x=[0: dx : 3∗ lambda ] ; % x axis [ m ]11 y=[0: dy : 3∗ lambda ] ; % y axis [ m ]1213 Nt=length ( t ) ; % # of t points [#]14 Nx=length ( x ) ; % # of x points [#]15 Ny=length ( y ) ; % # of y points [#]1617 E=zeros ( Nx , Ny , Nt ) ; % Electric field [ V/m ]1819 figure (1 )20 plot_handle=surf (x , y , E ( : , : , 1 ) ) ;21 shading flat22 xlabel ( ’ x [m] ’ )23 ylabel ( ’ y [m] ’ )24 zlabel ( ’E(x , n) ’ )25 h_title=title ( ’ ’ ) ;26 start_time=tic ; % start timer2728 f o r n=1:Nt29 f o r i=1:Nx30 f o r j=1:Ny31 E (i , j , n )=ones ( length ( y ( j ) ) , 1 ) ∗sin ( k∗x ( i )−2∗pi∗f∗t ( n ) ) ;32 end33 end34 % update plot35 set ( plot_handle , ’ZDATA’ , E ( : , : , n ) ) ;36 set ( h_title , ’ S t r ing ’ , [ ’ time= ’ num2str ( toc ( start_time ) ) ’ s ’ ] ) ;37 drawnow38 end1 %% Vectorized loop ( faster )2 f o r n=1:Nt3 i=1:Nx ;4 j=1:Ny ;5 E (i , j , n )=ones ( length ( y ( j ) ) , 1 ) ∗sin ( k∗x ( i )−2∗pi∗f∗t ( n ) ) ;6 . . .VectorizedFigure 5.2: An example to demonstrate how vectorizing for loops improves the speed ofmatlab code. A surface plot of a simple plane wave is sped up by replacing loops throughindividual values such as for i=1:Nx, to simply the whole vector i=1:Nx;.79Bibliography[1] M. Bethune-Waddell and K. J. Chau. “Simulations of radiation pressure experimentsnarrow down the energy and momentum of light in matter”. In: Reports on Progressin Physics 78.12 (2015), p. 122401.[2] J. C. Maxwell. A Treatise on Electricity and Magnetism. Vol. 2. Dover, 1891.[3] P. N. Lebedev. “Untersuchungen u¨ber die Druckkra¨fte des Lichtes[InvestigationsOver the compressive forces of light]”. In: Ann. Phys. 6 (1900), 433–458.[4] E. F. Nichols and G. F. Hull. “A Preliminary communication on the pressure of heatand light radiation”. In: Phys. Rev. (Series 1) 13 (1901), 307–320.[5] E. F. Nichols and G. F. Hull. “The pressure due to radiation. (Second Paper.)” In:Phys. Rev. (Series 1) 17 (1903), 26–50.[6] H. Minkowski. “Die Grundgleichungen fu¨r die elektromagnetischen Vorga¨nge in be-wegten Ko¨rpern”. In: Nachr. Ges. Wiss. Go¨ttingen (1908), 53–111.[7] M. Abraham. “Zur Elektrodynamik bewegter Ko¨rper[The electrodynamics of movingbodies]”. In: R. C. Circ. Mat. Palermo 28 (1909), 1–28.[8] A. Einstein and J. Laub. “U¨ber die elektromagnetischen Grundgleichungen fu¨r be-wegte Ko¨rper[About the electromagnetic fundamental equations for moving bodies]”.In: Ann. Phys. 26 (1908). Trans. A. Beck The Collected Papers of Albert Einstein,329–348, Princeton University Press, Princeton, 1989., 541–550.[9] R. M. Fano, L. J. Chu, and R. B. Adler. Electromagnetic fields, energy, and forces.New York: John Wiley & Sons, 1960.[10] L. J. Chu, H. Haus, and P. Penfield Jr. “The force density in polarizable and mag-netizable fluids”. In: Proc. IEEE 54 (1966), 920–935.[11] P. Penfield and H. A. Haus. Electrodynamics of Moving Media. Cambridge: M.I.T.Press, 1967.[12] Y. N. Obukhov and F. W. Hehl. “Electromagnetic energy-momentum and forces inmatter”. In: Phys. Lett. A 311 (2003), 277–284.80[13] R. V. Jones and J. C. S. Richards. “The Pressure of Radiation in a RefractingMedium”. In: Proc. R. Soc. Lond. 221 (1954), 480–498.[14] R. V. Jones and B. Leslie. “The Measurement of Optical Radiation Pressure inDispersive Media”. In: Proc. R. Soc. Lond. A 360 (1978), 347–363.[15] A. Ashkin and J. M. Dziedzic. “Radiation pressure on a free liquid surface”. In:Phys. Rev. Lett. 30 (1973), 139–142.[16] A. Casner and J-P. Delville. “Giant deformations of a liquid-liquid interface inducedby the optical radiation pressure”. In: Phys. Rev. Lett. 87 (2001), 054503.[17] A. Casner, J-P. Delville, and I. Brevik. “Asymmetric optical radiation pressure effectson liquid interfaces under intense illumination”. In: J. Opt. Soc. Am. B 20 (2003),2355–2362.[18] A. Casner and J-P. Delville. “Laser-Induced hydrodynamic instability of fluid inter-faces”. In: Phys. Rev. Lett. 90 (2003), 144503.[19] R. Wunenburger, A. Casner, and J-P. Delville. “Light-induced deformation and in-stability of a liquid interface. I. Statics”. In: Phys. Rev. E 73 (2006), 036314.[20] R. Wunenburger, A. Casner, and J-P. Delville. “Light-induced deformation and in-stability of a liquid interface. II. Dynamics”. In: Phys. Rev. E 73 (2006), 036315.[21] A. F. Gibson, M. F. Kimmitt, and A. C. Walker. “Photon drag in germanium”. In:Appl. Phys. Lett. 17 (1970), 75–77.[22] A. M. Danishevski˘i, A. A. Kastal’ski˘i, S. M. Ryvkin, and I. D. Yaroshetski˘i. “Drag-ging of free carriers by photons in direct interband transitions in semiconductors”.In: Sov. Phys. JETP 31 (1970), 292–295.[23] A. F. Gibson, M. F. Kimmitt, A. O. Koohian, D. E. Evans, and G. F. D. Levy. “Astudy of radiation pressure in a refractive medium by the photon drag effect”. In:Proc. R. Soc. Lond. A 370 (1980), 303–311.[24] G. K. Campbell, A. E. Leanhardt, J. Mun, M. Boyd, E. W. Streed, W. Ketterle,and D. E. Pritchard. “Photon recoil momentum in dispersive media”. In: Phys. Rev.Lett. 94 (2005), 170403.[25] F. N. H. Robinson. “Electromagnetic stress and momentum in matter”. In: Phys.Rep. 16 (1975), 313–354.[26] Z. Mikura. “Variational formulation of the electrodynamics of fluids and its applica-tion to the radiation pressure problem”. In: Phys. Rev. A 13 (1976), 2265–2275.[27] I. Brevik. “Experiments in phenomenological electrodynamics and the electromag-netic energy-momentum tensor”. In: Phys. Rep. 52 (1979), 133–201.81[28] B. A. Kemp and T. M. Grzegorczyk. “The observable pressure of light in dielectricfluids”. In: Opt. Lett. 36 (2011), 493–495.[29] N. L. Balazs. “The energy-momentum tensor of the electromagnetic field inside mat-ter”. In: Phys. Rev. 91 (1953), 408–411.[30] A. Einstein. “Das Prinzip von der Erhaltung der Schwerpunktsbewegung und dieTra¨gheit der Energie[The principle of conservation of the center of gravity and theinertia of energy]”. In: Ann. Phys. 20 (1906), 627–633.[31] L. D. Landau and E. M. Lifshitz. Electrodynamics of Continuous Media, SecondEdition. Vol. 8. Course of Theoretical Physics. 1984.[32] M. Kranys˘. “About the equivalence of Abraham’s and Minkowski’s electrodynam-ics”. In: Can. J. Phys. 57 (1979), 1022–1026.[33] R. N. C. Pfeifer, T. A. Nieminen, N. R. Heckenberg, and H. Rubinsztein-Dunlop.“Momentum of an electromagnetic wave in dielectric media”. In: Rev. Mod. Phys.79 (2007), pp. 1197–1216.[34] D. F. Nelson. “Momentum, pseudomomentum, and wave momentum: Toward re-solving the Minkowski-Abraham controversy”. In: Phys. Rev. A 44 (1991), 3985–3996.[35] E. A. Hinds and A. M. Barnett. “Momentum exchange between light and a singleatom: Abraham or Minkowski?” In: Phys. Rev. Lett. 102 (2009), 050403.[36] P. W. Milonni and R. W. Boyd. “Momentum of light in a dielectric medium”. In:Adv. Opt. Phot. 2 (2010), 519–553.[37] S. M. Barnett. “Resolution of the Abraham-Minkowski dilemma”. In: Phys. Rev.Lett. 104 (2010), 070401.[38] C. Baxter and R. Loudon. “Radiation pressure and the photon momentum in di-electrics”. In: J. Mod. Opt. 57 (2010), 830–842.[39] B. A. Kemp. “Resolution of the Abraham-Minkowski debate: Implications for theelectromagnetic wave theory of light in matter”. In: J. Appl. Phys. 109 (2011),111101.[40] R. V. Jones. “Radiation pressure of light in a dispersive medium”. In: Proc. R. Soc.Lond. A 360 (1978), 365–371.[41] R. Peierls. “The momentum of light in a refracting medium”. In: Proc. R. Soc. Lond.A 347 (1976), 475–491.[42] M. Mansuripur. “Radiation pressure and the linear momentum of the electromag-netic field”. In: Opt. Express 12 (2004), 5375–5401.82[43] M. Mansuripur. “Radiation pressure and the linear momentum of light in dispersivedielectric media”. In: Opt. Express 13 (2005), 2245–2250.[44] M. Mansuripur, A. R. Zakharian, and J. V. Moloney. “Radiation pressure on adielectric wedge”. In: Opt. Express 13 (2005), 2064–2074.[45] J. C. Garrison and R. Y. Chiao. “Canonical and kinetic forms of the electromagneticmomentum in an ad hoc quantization scheme for a dispersive dielectric”. In: Phys.Rev. A 70 (2004), 053826.[46] P. L. Saldanha. “Division of the momentum of electromagnetic waves in linear mediainto electromagnetic and material parts”. In: Opt. Express 18 (2011), 2258–2268.[47] I. Brevik and S. A˚. Ellingsen. “Electromagnetic momentum conservation in media”.In: Ann. Phys. 326 (2011), 754–769.[48] M. Mansuripur. “Momentum exchange effect”. In: Nature Photon. 7 (2013), 765–766.[49] M. V. Laue. “Zur Minkowskischen Elektrodynamik der bewegten Ko¨per.[The Minkowskielectrodynamics of moving bodies]”. In: Z. Physik 128 (1950), 387–394.[50] J. H. Poynting. “On the transfer of energy in the electromagnetic field”. In: Phil.Trans. R. Soc. 175 (1884), 343–361.[51] J. D. Jackson. Magnetic Actuators and Sensors, 2nd Edition. John Wiley and Sons,2014.[52] “Pull Force Case 3”. In: K&J Magnetics Inc. (2015). Accessed: 2015-11-01. url:https://www.kjmagnetics.com/largergraph.asp?CI=3\&pName=DAC.[53] N. L. Balazs. “On the propagation of energy in elastic media”. In: Proc. Phys. Soc.A 67 (1954), 726–727.[54] L. D. Landau and E. M. Lifshitz. The Classical Theory of Fields. Vol. 2. Course ofTheoretical Physics. 1951. isbn: 978-0-08-030275-1. doi: http://dx.doi.org/10.1016/B978-0-08-030275-1.50025-4.[55] J. D. Jackson. Classical Electrodynamics, Third Edition. John Wiley and Sons, 1999.[56] C. S. Helrich. The Classical Theory of Fields-Electromagnetism. Vol. 2. GraduateTexts in Physics. Springer, 2012.[57] M. Abraham. “Zur Frage der Symmetrie des elektroniagnetischen Spannungsten-sors[On the issue of symmetry of the electromagnetic stress tensor]”. In: Annal.Phys. 349 (1914), 537–544.[58] D. Griffths and I. M. Smith. Numerical methods for engineers. Chapman & Hall,2006.83[59] A. Taflove and S. C. Hagness. Computational Electrodynamics: The Finite-DifferenceTime-Domain Method, 3rd Ed. Norwood, MA: Artec House, 2005.[60] K. Yee. “Numerical solution of initial boundary value problems involving Maxwell’sequations in isotropic media”. In: vol. 13. 3. IEEE Transactions on Antennas andPropagation, 1966.[61] S. N. Makarov. Antenna and EM Modeling with Matlab. New York: Wiley & Sons,2002.[62] L. Yaxun and C. D. Sarris. “AMR-FDTD: a dynamically adaptive mesh refinementscheme for the finite-difference time-domain technique”. In: Antennas and Propaga-tion Society International Symposium, 2005 IEEE. Vol. 1A. 2005, 134–137.[63] G. Tryggvason. “A Front-tracking/Finite-Volume Navier-Stokes Solver for DirectNumerical Simulations of Multiphase Flows”. In: (2012). Unpublished. url: www3.nd.edu/~gtryggva/MultiphaseDNS/DNS-Solver.pd.[64] G. Tryggvason, B. Bunner, A. Esmaeeli, D. Juric, N. Al-Rawahi, W. Tauber W,J. Han, S. Nas, and J-Y. Jan. “A front-tracking method for the computations ofmultiphase flow”. In: J. Comp. Phys. 169 (2001), 708–759.[65] M. Mansuripur. Momentum in Classical Electrodynamics.: Bentham Books, 2011.[66] S. M. Barnett and R. Loudon. “The enigma of optical momentum in a medium”. In:Phil. Trans. R. Soc. A 368 (2010), 927–939.[67] K. J. Chau and H. Lezec. “Revisiting the Balazs thought experiment in the case ofa left-handed material: electromagnetic-pulse-induced displacement of a dispersive,dissipative negative-index slab”. In: Opt. Express 20 (2012), 10138–10162.[68] M. Sadiku. Elements of Electromagnetics. 4th ed. Oxford University Press, 2007.[69] D. M. Pozar. Microwave Engineering. 3rd ed. Wiley, 2005.[70] J. Blumm, A. Lindermann, M. Meyer, and C. Strasser. “Characterization of PTFEUsing Advanced Thermal Analysis Techniques”. In: Int. J. Themophys. 31 (2008),1919–1927.[71] I. Perepechko. Low-Temperature Properties of Polymers.: Elsevier, 2013.[72] M. Mansuripur. “Nature of electric and magnetic dipoles gleaned from the Poyntingtheorem and the Lorentz force law of classical electrodynamics”. In: Opt. Commun.284 (2011), 594–602.[73] A. Rohrbach and E. H. K. Stelzer. “Optical trapping of dielectric particles in arbi-trary fields”. In: J. Opt. Soc. Am. A 18 (2001), 839–853.84[74] M. Bell and S. E. Green. “On radiometer action and the pressure of radiation”. In:Proc. Phys. Soc. 45 (1933), 320–357.[75] M. Mansuripur. “Radiation pressure on submerged mirrors: implications for themomentum of light in dielectric media”. In: Opt. Express 15 (2007), 2677–2682.[76] M. Mansuripur. “Deducing radiation pressure on a submerged mirror from thedoppler shift”. In: Phys. Rev. A 85 (2012), 023807.[77] K. J. Webb. “Dependence of the radiation pressure on the background refractiveindex”. In: Phys. Rev. Lett. 111 (2013), 043602.[78] C. J. Sheppard and B. A. Kemp. “Optical pressure deduced from energy relationswithin relativistic formulations of electrodynamics”. In: Phys. Rev. A 89 (2014),013825.[79] E. Palik and G. Ghosh. Handbook on optical constants of solids. San Diego: AcademicPress, 1998.[80] P. B. Johnson and R. W. Christy. “Optical constants of the noble metals”. In: Phys.Rev. B 4 (1972), 4370–4379.[81] B. A. Kemp. “Macroscopic theory of optical momentum”. In: Prog. Opt. 60 (2015),437–488.[82] M. G. Burt and R. Peierls. “The momentum of a light Wave in a refracting medium”.In: Proc. R. Soc. Lond. A 333 (1973), 149–156.[83] N. G. C. Astrath, G. V. B. Lukasievicz, L. C. Malacarne and S. E. Bialkowski.“Surface deformation effects induced by radiation pressure and electrostriction forcesin dielectric solids”. In: Appl. Phys. Lett. 102 (2013), 231903.[84] J. P. Gordon. “Radiation Forces and Momenta in Dielectric”. In: Phys. Rev. A 8(1973), 14–21.[85] R. Loudon. “Radiation pressure and momentum in dielectrics”. In: Fortschr. Phys.52 (2004), 1134–1140.[86] N. G. C. Astrath, L. C. Malacarne, M. L. Baesso, G. V. B. Lukasievicz and S. E.Bialkowski. “Unravelling the effects of radiation forces in water”. In: Nat. Commun.5 (2014), 023826.[87] N. S. Aanensen, S. A˚. Ellingsen, and I. Brevik. “Theoretical considerations of laser-induced liquid-liquid interface deformation”. In: Phys. Scr. 87 (2013), 055402.[88] M. Mansuripur. “On the foundational equations of the classical theory of electrody-namics”. In: Resonance 18 (2013), 130–155.85[89] A. Casner. “De´formations, manipulations et instabilite´s d’interfaces liquides induitespar la pression de radiation d’une onde laser”. PhD thesis. University of Bordeaux,2002.[90] K. C. Neumann, E. H. Liou, G. F. Bergman, and S. M. Block. “Characterization ofphotodamage to Escherichia coli in optical traps”. In: Biophys J. 77 (1999), 2856–2863.[91] M. D. Wang, H. Yin, R. Landick, J. Gelles, and S. M. Block. “Stretching DNA withoptical tweezers”. In: Biophys. J. 72 (1997), 1335–1346.[92] P. H. Jones, O. M. Marago, and G. Volpe. Optical Tweezers: Principles and Appli-cations. Cambridge University Press, 2016.[93] D. G. Grier. “A revolution in optical manipulation”. In: Nature 424 (2003), 810–816.[94] J. P. Barton, D. R. Alexander and S. A. Schaub. “Theoretical determination of netradiation force and torque for a spherical particle illuminated by a focused laserbeam”. In: J. Appl. Phys. 66 (1989), 4594–4602.[95] R. N. C. Pfeifer, T. A. Nieminen, N. R. Heckenberg and H. Rubinsztein-Dunlop.“Optical tweezers and paradoxes in electromagnetism”. In: J. Opt. 13 (2011), 044017.[96] A. Ashkin. “Acceleration and trapping of particles by radiation pressure”. In: Phys.Rev. Lett. 24 (1970), 156–159.[97] P. L. Marston. “Axial radiation force of a Bessel beam on a sphere and directionreversal of the force”. In: J. Acoust. Soc. Am. 120 (2006), 3518.[98] S. Sukhov and A. Dogariu. “On the concept of “tractor beams””. In: Opt. Lett. 35(2010), 3847–3849.[99] S. Sukhov and A. Dogariu. “Negative nonconservative Forces: Optical “TractorBeams” for arbitrary objects”. In: Phys. Rev. Lett. 107 (2011), 203602.[100] A. Novitsky, C-W. Qiu and H. Wang. “Single gradientless light beam drags particlesas tractor beams”. In: Phys. Rev. Lett. 107 (2011), 203601.[101] J. Chen, Y. Wang, B. Jia, T. Geng, X. Li, L. Feng, W. Qian, B. Liang, X. Zhang,M. Gu, and S. Zhuang. “Observation of the inverse Doppler effect in negative-indexmaterials at optical frequencies”. In: Nat. Phot. 5 (2011), 239–245.[102] D. B. Ruffner and D. G. Grier. “Optical conveyors: a class of active tractor beams”.In: Phys. Rev. Lett. 109 (2012), 163903.[103] V. Kajorndejnukul, W. Ding, S. Sukhov, C-W. Qiu and A. Dogariu. “Linear momen-tum increase and negative optical forces at dielectric interface”. In: Nature Photonics7 (2013), 787–790.86[104] O. Brzobohaty´, V. Kara´sek, M. Sˇiler, L. Chva´tal, T. Cˇizˇma´r and P. Zema´nek. “Ex-perimental demonstration of optical transport, sorting and self-arrangement using a‘tractor beam’”. In: Nature Photon. 7 (2013), 123–127.[105] V. G. Veselago. “The Electrodynamics of Substances with Simultaneously Negativevalues of ε and µ”. In: Sov. Phys. Usp. 10 (1968), 509–514.[106] R. A. Shelby, D. R. Smith, and S. Schultz. “Experimental Verification of a NegativeIndex of Refraction”. In: Science 292 (5514 2001), 77–79.[107] S. M. Rytov. “Some Theorems on the group velocity of electromagnetic waves”. In:J. Exp. Th. Phys. 17 (1947), 930–936.[108] V. P. Makarov, and A. A. Rukhadze. “Negative group velocity electromagnetic wavesand the energy-momentum tensor”. In: 54 (12 2011), pp. 1285–1296.87
- Library Home /
- Search Collections /
- Open Collections /
- Browse Collections /
- UBC Theses and Dissertations /
- Simulations of radiation pressure experiments
Open Collections
UBC Theses and Dissertations
Featured Collection
UBC Theses and Dissertations
Simulations of radiation pressure experiments Bethune-Waddell, Maximilien 2016
pdf
Page Metadata
Item Metadata
Title | Simulations of radiation pressure experiments |
Creator |
Bethune-Waddell, Maximilien |
Publisher | University of British Columbia |
Date Issued | 2016 |
Description | Electromagnetic radiation, such as light, can exert a force on objects. This phenomena is referred to as radiation pressure. A popular example would be the radiation pressure from the sun which is presently the only force pushing NASA's Kepler II satellite out of our solar system. The most important material property relevant to radiation pressure is the refractive index n. The vacuum of space around the Kepler satellite has a refractive index of n=1. For this situation our ability to model radiation pressure is quite certain. There has been a debate over the last century on how to model radiation pressure in the presence of matter, where the refractive index of the surrounding medium can be n >1. Five prevalent models exist for electrodynamics, namely the Abraham, Minkowksi, Einstein-Laub, Chu, or the Amperian formulations. Various electrodynamic theories presented over the last century have differences that only arise in exotic situations such as objects in liquid environments, material moving at relativistic speeds, and short timescales difficult to measure. Because technology has now advanced to potentially include these situations, this thesis addresses the renewed attention this debate requires. Each of these models are tested against two criteria using a simulation tool that solves for both electromagnetic and fluid dynamic phenomena. First, their compliance with the conservation laws of energy, momentum, and center-of-mass velocity. Second, their accord with experimental results. A simulation environment is used to calculate the conserved quantities and observable effects that each model predicts under various experimental conditions. The simulation tool is a Matlab script that allows us to consider and compare the results of many experiments simultaneously. Five significant experiments are analyzed in this thesis: the radiation pressure observed on metallic and dielectric mirrors, the deformation of a water-air interface, the deformation of a fluid-fluid interface, the displacement of polystyrene beads submerged in water, and the displacement of an oil droplet on a water surface. This thesis reaches the conclusion that not enough data is available from past experiments to verify a single electrodynamic theory. Our work suggests simple new experiments that could. |
Genre |
Thesis/Dissertation |
Type |
Text |
Language | eng |
Date Available | 2016-09-01 |
Provider | Vancouver : University of British Columbia Library |
Rights | Attribution-NoDerivatives 4.0 International |
IsShownAt | 10.14288/1.0313408 |
URI | http://hdl.handle.net/2429/59024 |
Degree |
Master of Applied Science - MASc |
Program |
Engineering |
Affiliation |
Applied Science, Faculty of Engineering, School of (Okanagan) |
Degree Grantor | University of British Columbia |
GraduationDate | 2016-11 |
Campus |
UBCO |
Scholarly Level | Graduate |
Rights URI | http://creativecommons.org/licenses/by-nd/4.0/ |
AggregatedSourceRepository | DSpace |
Download
- Media
- 24-ubc_2016_november_Bethune-Waddell_Maximilien.pdf [ 5.51MB ]
- Metadata
- JSON: 24-1.0313408.json
- JSON-LD: 24-1.0313408-ld.json
- RDF/XML (Pretty): 24-1.0313408-rdf.xml
- RDF/JSON: 24-1.0313408-rdf.json
- Turtle: 24-1.0313408-turtle.txt
- N-Triples: 24-1.0313408-rdf-ntriples.txt
- Original Record: 24-1.0313408-source.json
- Full Text
- 24-1.0313408-fulltext.txt
- Citation
- 24-1.0313408.ris
Full Text
Cite
Citation Scheme:
Usage Statistics
Share
Embed
Customize your widget with the following options, then copy and paste the code below into the HTML
of your page to embed this item in your website.
<div id="ubcOpenCollectionsWidgetDisplay">
<script id="ubcOpenCollectionsWidget"
src="{[{embed.src}]}"
data-item="{[{embed.item}]}"
data-collection="{[{embed.collection}]}"
data-metadata="{[{embed.showMetadata}]}"
data-width="{[{embed.width}]}"
async >
</script>
</div>
Our image viewer uses the IIIF 2.0 standard.
To load this item in other compatible viewers, use this url:
http://iiif.library.ubc.ca/presentation/dsp.24.1-0313408/manifest