Maxwell-Random Debye System Distributions of Relaxation Times Distributions of Relaxation Times As early as 1897, experiments by Drude demonstrated that some materials exhibited “anomalous dispersion, i.e, the decrease in the dielectric constant with increasing frequency” [Debye 1929]. In 1907, von Schweidler observed the need for multiple relaxation times. Around the same time, Debye’s papers appeared (in German) defining and quantifying the relaxation time. Analogous to the Maxwell-Wiechert model of viscoelasticity from 1893. In 1913, Wagner proposed a continuous distribution of relaxation times. In 1927, Debye invited to US (and translated his works into English). In 1927, K.S. Cole studied electrical properties of biological systems (during his postdoc at Harvard). Was invited to visit Debye in Germany in 1928. N. L. Gibson (Oregon State) Maxwell-PC Dispersive ICERM 2018 12 / 72
Maxwell-Random Debye System Distributions of Relaxation Times In 1940, R.H. Cole began PhD at Harvard and collaborated with brother on a graphical method to test the Debye model, as well as a heuristic fix. N. L. Gibson (Oregon State) Maxwell-PC Dispersive ICERM 2018 13 / 72
Maxwell-Random Debye System Distributions of Relaxation Times In 1940, R.H. Cole began PhD at Harvard and collaborated with brother on a graphical method to test the Debye model, as well as a heuristic fix. In 1950, D.W. Davidson and R.H. Cole discovered materials that are not well-represented by the Cole-Cole model, proposed a different heuristic. N. L. Gibson (Oregon State) Maxwell-PC Dispersive ICERM 2018 13 / 72
Maxwell-Random Debye System Distributions of Relaxation Times In 1940, R.H. Cole began PhD at Harvard and collaborated with brother on a graphical method to test the Debye model, as well as a heuristic fix. In 1950, D.W. Davidson and R.H. Cole discovered materials that are not well-represented by the Cole-Cole model, proposed a different heuristic. In 1967, Havriliak and Negami combined the two heuristics to form a generalized model. N. L. Gibson (Oregon State) Maxwell-PC Dispersive ICERM 2018 13 / 72
Maxwell-Random Debye System Distributions of Relaxation Times In 1940, R.H. Cole began PhD at Harvard and collaborated with brother on a graphical method to test the Debye model, as well as a heuristic fix. In 1950, D.W. Davidson and R.H. Cole discovered materials that are not well-represented by the Cole-Cole model, proposed a different heuristic. In 1967, Havriliak and Negami combined the two heuristics to form a generalized model. One can show that the Cole-Cole model (and extensions) corresponds to a continuous distribution of relaxation times “... it is possible to calculate the necessary distribution function by the method of Fuoss and Kirkwood.” [Cole-Cole1941]. N. L. Gibson (Oregon State) Maxwell-PC Dispersive ICERM 2018 13 / 72
Maxwell-Random Debye System Distributions of Relaxation Times In 1940, R.H. Cole began PhD at Harvard and collaborated with brother on a graphical method to test the Debye model, as well as a heuristic fix. In 1950, D.W. Davidson and R.H. Cole discovered materials that are not well-represented by the Cole-Cole model, proposed a different heuristic. In 1967, Havriliak and Negami combined the two heuristics to form a generalized model. One can show that the Cole-Cole model (and extensions) corresponds to a continuous distribution of relaxation times “... it is possible to calculate the necessary distribution function by the method of Fuoss and Kirkwood.” [Cole-Cole1941]. N. L. Gibson (Oregon State) Maxwell-PC Dispersive ICERM 2018 13 / 72
Maxwell-Random Debye System Distributions of Relaxation Times In 1940, R.H. Cole began PhD at Harvard and collaborated with brother on a graphical method to test the Debye model, as well as a heuristic fix. In 1950, D.W. Davidson and R.H. Cole discovered materials that are not well-represented by the Cole-Cole model, proposed a different heuristic. In 1967, Havriliak and Negami combined the two heuristics to form a generalized model. One can show that the Cole-Cole model (and extensions) corresponds to a continuous distribution of relaxation times “... it is possible to calculate the necessary distribution function by the method of Fuoss and Kirkwood.” [Cole-Cole1941]. y 2 α + 2 y α cos( πα ) + 1 F ( y ) = y αβ � � sin( βθ ) /π − β/ 2 , where y = τ/τ 0 and θ is defined implictly by ( y α + cos ( πα )) tan ( θ ) = sin( πα ) . N. L. Gibson (Oregon State) Maxwell-PC Dispersive ICERM 2018 13 / 72
Maxwell-Random Debye System Distributions of Relaxation Times [Garrappa2016] N. L. Gibson (Oregon State) Maxwell-PC Dispersive ICERM 2018 14 / 72
Maxwell-Random Debye System Distributions of Relaxation Times Figure: Relaxation Time Distribution for CC model [Garrappa2016]. N. L. Gibson (Oregon State) Maxwell-PC Dispersive ICERM 2018 15 / 72
Maxwell-Random Debye System Fit to dry skin data with uniform distribution True Data 3 10 Debye (27.79) Cole−Cole (10.4) Uniform (13.60) ε 2 10 2 4 6 8 10 10 10 10 10 10 f (Hz) Figure: Real part of ǫ ( ω ), ǫ , or the permittivity [REU2008]. N. L. Gibson (Oregon State) Maxwell-PC Dispersive ICERM 2018 16 / 72
Maxwell-Random Debye System Fit to dry skin data with uniform distribution 1 True Data 10 Debye (27.79) Cole−Cole (10.4) Uniform (13.60) 0 10 −1 10 σ −2 10 −3 10 2 4 6 8 10 10 10 10 10 10 f (Hz) Figure: Imaginary part of ǫ ( ω ) /ω , σ , or the conductivity [REU2008]. N. L. Gibson (Oregon State) Maxwell-PC Dispersive ICERM 2018 17 / 72
Maxwell-Random Debye System Distribution of Relaxation Times Distributions of Parameters To account for the effect of distributions of parameters q , consider the following polydispersive DRF � h ( t , x ; F ) = g ( t , x ; q ) dF ( q ) , Q where Q is some admissible set and F ∈ P ( Q ). Then the polarization becomes: � t P ( t , x ; F ) = h ( t − s , x ; F ) E ( s , x ) ds . 0 Alternatively we can define the random polarization P ( t , x ; q ) to satisfy � P ( t , x ; q ) dF ( q ) . P ( t , x ; F ) = Q N. L. Gibson (Oregon State) Maxwell-PC Dispersive ICERM 2018 18 / 72
Maxwell-Random Debye System Distribution of Relaxation Times Random Polarization In the case of relaxation polarization, the random polarization P ( t , x ; τ ) solves τ ˙ P + P = ǫ 0 ǫ d E where τ is a random variable with PDF f ( τ ), for example, 1 f ( τ ) = τ b − τ a for a uniform distribution. N. L. Gibson (Oregon State) Maxwell-PC Dispersive ICERM 2018 19 / 72
Maxwell-Random Debye System Distribution of Relaxation Times Random Polarization In the case of relaxation polarization, the random polarization P ( t , x ; τ ) solves τ ˙ P + P = ǫ 0 ǫ d E where τ is a random variable with PDF f ( τ ), for example, 1 f ( τ ) = τ b − τ a for a uniform distribution. The electric field depends on the macroscopic polarization, which we take to be the expected value of the random polarization at each point ( t , x ) � τ b P ( t , x ; τ ) f ( τ ) d τ. P ( t , x ; F ) = τ a N. L. Gibson (Oregon State) Maxwell-PC Dispersive ICERM 2018 19 / 72
Maxwell-Random Debye System Polynomial Chaos Polynomial Chaos Apply Polynomial Chaos (PC) method to approximate each spatial component of the random polarization τ ˙ P + P = ǫ 0 ǫ d E , τ = τ ( ξ ) = τ r ξ + τ m , ξ ∼ F (with ξ mean 0 and variance 1) resulting in ( τ r M + τ m I ) ˙ α + � � α = ǫ 0 ǫ d E ˆ e 1 or A ˙ α = � � α + � f . N. L. Gibson (Oregon State) Maxwell-PC Dispersive ICERM 2018 20 / 72
Maxwell-Random Debye System Polynomial Chaos Polynomial Chaos Apply Polynomial Chaos (PC) method to approximate each spatial component of the random polarization τ ˙ P + P = ǫ 0 ǫ d E , τ = τ ( ξ ) = τ r ξ + τ m , ξ ∼ F (with ξ mean 0 and variance 1) resulting in ( τ r M + τ m I ) ˙ � α + � α = ǫ 0 ǫ d E ˆ e 1 or A ˙ α = � � α + � f . The electric field depends on the macroscopic polarization, the expected value of the random polarization at each point ( t , x ), which is P ( t , x ; F ) = E [ P ] ≈ α 0 ( t , x ) . Note that A is positive definite if τ r < τ m since λ ( M ) ∈ ( − 1 , 1). N. L. Gibson (Oregon State) Maxwell-PC Dispersive ICERM 2018 20 / 72
Maxwell-Random Debye System Polynomial Chaos Maximum Difference Calculated for different values of p and r 2 10 r = 1 . 00 τ 0 10 r = 0 . 75 τ r = 0 . 50 τ r = 0 . 25 τ −2 10 Maximum Error −4 10 −6 10 −8 10 −10 10 −12 10 0 1 2 3 4 5 6 7 8 9 10 p N. L. Gibson (Oregon State) Maxwell-PC Dispersive ICERM 2018 21 / 72
Maxwell-Random Debye System Inverse Problem Numerical Results 11 Distributions, noise = 0.1, refinement = 1, perturb = −0.8 x 10 6 Initial J=983.713 Optimal J=1.25869 Actual J=1.25879 5 4 3 2 1 0 −12 −11 −10 10 10 10 Comparison of initial to final distribution [Armentrout-G., 2011]. N. L. Gibson (Oregon State) Maxwell-PC Dispersive ICERM 2018 22 / 72
Maxwell-Random Debye System Inverse Problem Numerical Results Comparison, noise = 0.1, refinement = 1, perturb = −0.8 100 Data Initial J=983.713 Optimal J=1.25869 80 Actual J=1.25879 60 40 20 E 0 −20 −40 −60 −80 0 0.2 0.4 0.6 0.8 1 1.2 1.4 1.6 t −9 x 10 Comparison of simulations to data [Armentrout-G., 2011]. N. L. Gibson (Oregon State) Maxwell-PC Dispersive ICERM 2018 23 / 72
Maxwell-Random Debye System Inverse Problem Numerical Results 3.5 × 10 11 Distributions, noise = 0.1, refinement = 1,test = 5 Initial J=9.59699 Optimal J=3.02332 Actual J=3.02348 3 2.5 2 1.5 1 0.5 0 0 0.2 0.4 0.6 0.8 1 1.2 × 10 -11 N. L. Gibson (Oregon State) Maxwell-PC Dispersive ICERM 2018 24 / 72
Stability and Dispersion Analyses Debye Dispersion Relation (Deterministic) Maxwell-Debye System Combining Maxwell’s Equations, Constitutive Laws, and the Debye model, we have ∂ H µ 0 ∂ t = −∇ × E , (1a) ∂ E ∂ t = ∇ × H − ǫ 0 ǫ d E + 1 ǫ 0 ǫ ∞ τ P − J , (1b) τ τ ∂ P ∂ t = ǫ 0 ǫ d E − P . (1c) N. L. Gibson (Oregon State) Maxwell-PC Dispersive ICERM 2018 25 / 72
Stability and Dispersion Analyses Debye Dispersion Relation Assuming a solution to (1) of the form E = E 0 exp ( i ( ω t − k · x )), the following relation must hold. Debye Dispersion Relation The dispersion relation for the Maxwell-Debye system is given by ω 2 c 2 ǫ ( ω ) = � k � 2 where the complex permittivity is given by � 1 � ǫ ( ω ) = ǫ ∞ + ǫ d 1 + i ωτ Here, k is the wave vector and c = 1 / √ µ 0 ǫ 0 is the speed of light. N. L. Gibson (Oregon State) Maxwell-PC Dispersive ICERM 2018 26 / 72
Stability and Dispersion Analyses Random Debye Dispersion Relation Maxwell-Random Debye system In a polydispersive Debye material, we have ∂ H µ 0 ∂ t = −∇ × E , (2a) ∂ E ∂ t = ∇ × H − ∂ P ǫ 0 ǫ ∞ ∂ t − J (2b) τ ∂ P ∂ t + P = ǫ 0 ǫ d E (2c) with � τ b P ( t , x ; τ ) dF ( τ ) . P ( t , x ; F ) = τ a N. L. Gibson (Oregon State) Maxwell-PC Dispersive ICERM 2018 27 / 72
Stability and Dispersion Analyses Random Debye Dispersion Relation Theorem (G., 2015) The dispersion relation for the system (14) is given by ω 2 c 2 ǫ ( ω ) = � k � 2 where the expected complex permittivity is given by � 1 � ǫ ( ω ) = ǫ ∞ + ǫ d E . 1 + i ωτ Where k is the wave vector and c = 1 / √ µ 0 ǫ 0 is the speed of light. N. L. Gibson (Oregon State) Maxwell-PC Dispersive ICERM 2018 28 / 72
Stability and Dispersion Analyses Random Debye Dispersion Relation Theorem (G., 2015) The dispersion relation for the system (14) is given by ω 2 c 2 ǫ ( ω ) = � k � 2 where the expected complex permittivity is given by � 1 � ǫ ( ω ) = ǫ ∞ + ǫ d E . 1 + i ωτ Where k is the wave vector and c = 1 / √ µ 0 ǫ 0 is the speed of light. Note: for a uniform distribution on [ τ a , τ b ], this has an analytic form since 1 + ( ωτ ) 2 �� τ = τ a � � � 1 1 arctan( ωτ ) + i 1 � = 2 ln . E 1 + i ωτ ω ( τ b − τ a ) τ = τ b The exact dispersion relation can be compared with a discrete dispersion relation to determine the amount of dispersion error. N. L. Gibson (Oregon State) Maxwell-PC Dispersive ICERM 2018 28 / 72
Stability and Dispersion Analyses Debye Stability 2D Maxwell-Debye Transverse Electric (TE) curl equations For simplicity in exposition and to facilitate analysis, we reduce the Maxwell-Debye model to two spatial dimensions (we make the assumption that fields do not exhibit variation in the z direction). ∂ H µ 0 ∂ t = − curl E , (3a) ∂ E ∂ t = curl H − ǫ 0 ǫ d E + 1 τ P − J , ǫ 0 ǫ ∞ (3b) τ τ ∂ P ∂ t = ǫ 0 ǫ d E − P , (3c) where E = ( E x , E y ) T , P = ( P x , P y ) T and H z = H . � T � Note curl U = ∂ U y ∂ x − ∂ U x ∂ V ∂ y , − ∂ V ∂ y and curl V = . ∂ x N. L. Gibson (Oregon State) Maxwell-PC Dispersive ICERM 2018 29 / 72
Stability and Dispersion Analyses Debye Stability Stability Estimates for Maxwell-Debye System is well-posed since solutions satisfy the following stability estimate. Theorem (Li2010) Let D ⊂ R 2 , and let H, E , and P be the solutions to (the weak form of) the 2D Maxwell-Debye TE system with PEC boundary conditions. Then the system exhibits energy decay E ( t ) ≤ E (0) , ∀ t ≥ 0 where the energy is defined by 2 E ( t ) 2 = �√ µ 0 H ( t ) � 2 2 + �√ ǫ 0 ǫ ∞ E ( t ) � 2 � 1 � � � 2 + √ ǫ 0 ǫ d P ( t ) � � � � 2 and � · � 2 is the L 2 ( D ) norm. N. L. Gibson (Oregon State) Maxwell-PC Dispersive ICERM 2018 30 / 72
Stability and Dispersion Analyses Random Debye Stability We introduce the random Hilbert space V F = ( L 2 (Ω) ⊗ L 2 ( D )) 2 equipped with an inner product and norm as follows ( u , v ) F = E [( u , v ) 2 ] , � u � 2 F = E [ � u � 2 2 ] . The weak formulation of the 2D Maxwell-Random Debye TE system is � ∂ H � � − 1 � ∂ t , v = curl E , v , (4) µ 0 2 2 � ∂ E � � ∂ P � = ( H , curl u ) 2 − ǫ 0 ǫ ∞ ∂ t , u ∂ t , u , (5) 2 2 � ∂ P � � ǫ 0 ǫ d � 1 � � F − τ P , w ∂ t , w = E , w , (6) τ F F for v ∈ L 2 ( D ), u ∈ H 0 ( curl , D ) 2 , and w ∈ V F . N. L. Gibson (Oregon State) Maxwell-PC Dispersive ICERM 2018 31 / 72
Stability and Dispersion Analyses Random Debye Stability Stability Estimates for Maxwell-Random Debye System is well-posed since solutions satisfy the following stability estimate. Theorem (G., 2015) Let D ⊂ R 2 , and let H, E , and P be the solutions to the weak form of the 2D Maxwell-Random Debye TE system with PEC boundary conditions. Then the system exhibits energy decay E ( t ) ≤ E (0) , ∀ t ≥ 0 where the energy is defined by 2 2 + �√ ǫ 0 ǫ ∞ E ( t ) � 2 � � E ( t ) 2 = �√ µ 0 H ( t ) � 2 1 � � 2 + P ( t ) . √ ǫ 0 ǫ d � � � � F N. L. Gibson (Oregon State) Maxwell-PC Dispersive ICERM 2018 32 / 72
Stability and Dispersion Analyses Random Debye Stability Proof: (for 2D) By choosing v = H , u = E , and w = P in the weak form, and adding all three equations into the time derivative of the definition of E 2 , we obtain d E 2 ( t ) 1 � ǫ 0 ǫ d � 1 � � � � � � = − 2 − τ P , E curl E , H 2 + H , curl E E , E F + 2 dt τ F � 1 1 � � � τ E , P F − ǫ 0 ǫ d τ P , P + F � 1 � 1 1 � 1 � � � = − ǫ 0 ǫ d τ E , E F + 2 τ P , E F − τ P , P ǫ 0 ǫ d F 2 � � = − 1 1 � � τ ( P − ǫ 0 ǫ d E ) . � � ǫ 0 ǫ d � � F 2 � � d E ( t ) − 1 1 � � τ ( P − ǫ 0 ǫ d E ) ≤ 0 . = � � dt ǫ 0 ǫ d E ( t ) � � F N. L. Gibson (Oregon State) Maxwell-PC Dispersive ICERM 2018 33 / 72
Stability and Dispersion Analyses Maxwell-PC Debye Stability Maxwell-PC Debye Replace the Debye model with the PC approximation. In two dimensions we have the 2D Maxwell-PC Debye TE scalar equations ∂ H ∂ t = ∂ E x ∂ y − ∂ E y µ 0 ∂ x , (7a) ∂ E x ∂ t = ∂ H ∂ y − ∂α 0 , x ǫ 0 ǫ ∞ , (7b) ∂ t ∂ E y ∂ x − ∂α 0 , y ∂ t = − ∂ H ǫ 0 ǫ ∞ , (7c) ∂ t A ˙ α x = � � α x + � f x , (7d) A ˙ α y = � α y + � � f y . (7e) where � e 1 and � α y ] T . f x = ǫ 0 ǫ d E x ˆ f y = ǫ 0 ǫ d E y ˆ e 1 . Denote � α = [ � α x , � N. L. Gibson (Oregon State) Maxwell-PC Dispersive ICERM 2018 34 / 72
Stability and Dispersion Analyses Maxwell-PC FDTD Finite Difference Time Domain (FDTD) We now choose a discretization of the Maxwell-PC Debye model. Note that any scheme can be used independent of the spectral approach in random space employed here. The Yee Scheme (FDTD) This gives an explicit second order accurate scheme in time and space. It is conditionally stable with the CFL condition ν := c ∆ t 1 ≤ √ h d where ν is called the Courant number and c ∞ = 1 / √ µ 0 ǫ 0 ǫ ∞ is the fastest wave speed and d is the spatial dimension, and h is the (uniform) spatial step. The Yee scheme can exhibit numerical dispersion and dissipation. N. L. Gibson (Oregon State) Maxwell-PC Dispersive ICERM 2018 35 / 72
Stability and Dispersion Analyses Maxwell-PC FDTD Discrete Debye Dispersion Relation (Petropolous1994) showed that for the Yee scheme applied to the Maxwell-Debye, the discrete dispersion relation can be written ω 2 c 2 ǫ ∆ ( ω ) = K 2 ∆ ∆ where the discrete complex permittivity is given by � 1 � ǫ ∆ ( ω ) = ǫ ∞ + ǫ d 1 + i ω ∆ τ ∆ with discrete (mis-)representations of ω and τ given by ω ∆ = sin ( ω ∆ t / 2) , τ ∆ = sec( ω ∆ t / 2) τ . ∆ t / 2 N. L. Gibson (Oregon State) Maxwell-PC Dispersive ICERM 2018 36 / 72
Stability and Dispersion Analyses Maxwell-PC FDTD Discrete Debye Dispersion Relation (cont.) The quantity K ∆ is given by K ∆ = sin ( k ∆ z / 2) ∆ z / 2 in 1D and is related to the symbol of the discrete first order spatial difference operator by iK ∆ = F ( D 1 , ∆ z ) . In this way, we see that the left hand side of the discrete dispersion relation ω 2 c 2 ǫ ∆ ( ω ) = K 2 ∆ ∆ is unchanged when one moves to higher order spatial derivative approximations [Bokil-G,2012] or even higher spatial dimension [Bokil-G,2013]. N. L. Gibson (Oregon State) Maxwell-PC Dispersive ICERM 2018 37 / 72
Stability and Dispersion Analyses Maxwell-PC FDTD h , τ E y Let τ E x h , τ H h be the sets of spatial grid points on which the E x , E y , and H fields, respectively, will be discretized. The discrete L 2 grid norms are defined as L − 1 J − 1 � 2 , j | 2 + | V y ℓ, j + 1 2 | 2 � � V � 2 � � E = ∆ x ∆ y | V x ℓ + 1 , (8) ℓ =0 j =0 L − 1 J − 1 � U � 2 � � 2 | 2 , H = ∆ x ∆ y | U ℓ + 1 (9) 2 , j + 1 ℓ =0 j =0 with corresponding inner products. Each component α k is discretized on with discrete L 2 grid norm × τ E y τ E x h h p α � 2 � � α k � 2 � � α = E , k =0 with a corresponding inner product p � � α , � � ( � β ) α = α k , β k E . k =0 N. L. Gibson (Oregon State) Maxwell-PC Dispersive ICERM 2018 38 / 72
Stability and Dispersion Analyses PC-Debye FDTD Stability Energy Decay and Stability Energy decay implies that the method is stable and hence convergent. Theorem (G., 2015) For n ≥ 0 , let U n = [ H n − 1 0 , y , . . . ] T be the solutions of 2 , E n x , E n y , α n 0 , x , . . . , α n the 2D Maxwell-PC Debye TE FDTD scheme with PEC boundary conditions. If the usual CFL condition for Yee scheme is satisfied √ c ∞ ∆ t ≤ h / 2 , then there exists the energy decay property E n +1 ≤ E n h h where the discrete energy is given by 2 � √ µ 0 H H + ||√ ǫ 0 ǫ ∞ E n || 2 � � � � 2 1 � � n � � h ) 2 = ( E n � � α n � � E + √ ǫ 0 ǫ d � . � � � � � � � � � � � � � � � α N. L. Gibson (Oregon State) Maxwell-PC Dispersive ICERM 2018 39 / 72
Stability and Dispersion Analyses PC-Debye FDTD Stability Energy Decay and Stability Energy decay implies that the method is stable and hence convergent. Theorem (G., 2015) For n ≥ 0 , let U n = [ H n − 1 0 , y , . . . ] T be the solutions of 2 , E n x , E n y , α n 0 , x , . . . , α n the 2D Maxwell-PC Debye TE FDTD scheme with PEC boundary conditions. If the usual CFL condition for Yee scheme is satisfied √ c ∞ ∆ t ≤ h / 2 , then there exists the energy decay property E n +1 ≤ E n h h where the discrete energy is given by 2 � √ µ 0 H H + ||√ ǫ 0 ǫ ∞ E n || 2 � � � � 2 1 � � n � � h ) 2 = ( E n � � α n � � E + √ ǫ 0 ǫ d � . � � � � � � � � � � � � � � � α 2 ] = � E [ P ] 2 + Var ( P ) � 2 Note: �P� 2 F = E [ �P� 2 α � 2 2 ≈ � � α . N. L. Gibson (Oregon State) Maxwell-PC Dispersive ICERM 2018 39 / 72
Stability and Dispersion Analyses PC-Debye FDTD Stability Energy Decay and Stability (cont.) Proof. First, showing that this is a discrete energy, i.e., a positive definite function of the solution, involves recognizing that 1 h ) 2 = µ 0 � H α n − E ˆ α n − E ˆ n � 2 ( E n H + ǫ 0 ǫ ∞ ( E n , A h E n ) E + e 1 , A − 1 ( � ( � e 1 )) α ǫ 0 ǫ d with A h positive definite when the CFL condition is satisfied, and A − 1 is always positive definite (eigenvalues between τ m − τ r and τ m + τ r ). N. L. Gibson (Oregon State) Maxwell-PC Dispersive ICERM 2018 40 / 72
Stability and Dispersion Analyses PC-Debye FDTD Stability Energy Decay and Stability (cont.) Proof. First, showing that this is a discrete energy, i.e., a positive definite function of the solution, involves recognizing that 1 h ) 2 = µ 0 � H α n − E ˆ α n − E ˆ n � 2 ( E n H + ǫ 0 ǫ ∞ ( E n , A h E n ) E + e 1 , A − 1 ( � ( � e 1 )) α ǫ 0 ǫ d with A h positive definite when the CFL condition is satisfied, and A − 1 is always positive definite (eigenvalues between τ m − τ r and τ m + τ r ). The rest follows the proof for the deterministic case [Bokil-G, 2014] to show � � E n +1 2 − E n � � 2 1 n + 1 n + 1 h h � � = − � ǫ 0 ǫ d E 2 ˆ e 1 − � 2 A − 1 . (10) α � � E n +1 ∆ t ǫ 0 ǫ d + E n � h h N. L. Gibson (Oregon State) Maxwell-PC Dispersive ICERM 2018 40 / 72
Stability and Dispersion Analyses PC-Debye FDTD Dispersion Relation Theorem (G., 2015) The discrete dispersion relation for the Maxwell-PC Debye FDTD scheme is given by ω 2 c 2 ǫ ∆ ( ω ) = K 2 ∆ ∆ where the discrete expected complex permittivity is given by 1 ( I + i ω ∆ A ∆ ) − 1 ˆ e T ǫ ∆ ( ω ) := ǫ ∞ + ǫ d ˆ e 1 and the discrete PC matrix is given by A ∆ := sec( ω ∆ t / 2) A . The definitions of the parameters ω ∆ and K ∆ are the same as before. Recall the exact complex permittivity is given by � 1 � ǫ ( ω ) = ǫ ∞ + ǫ d E 1 + i ωτ . N. L. Gibson (Oregon State) Maxwell-PC Dispersive ICERM 2018 41 / 72
Stability and Dispersion Analyses PC-Debye FDTD Dispersion Analysis Dispersion Error We define the phase error Φ for a scheme applied to a model to be � � k EX − k ∆ � � Φ = � , (11) � � k EX � where the numerical wave number k ∆ is implicitly determined by the corresponding dispersion relation and k EX is the exact wave number for the given model. N. L. Gibson (Oregon State) Maxwell-PC Dispersive ICERM 2018 42 / 72
Stability and Dispersion Analyses PC-Debye FDTD Dispersion Analysis Dispersion Error We define the phase error Φ for a scheme applied to a model to be � � k EX − k ∆ � � Φ = � , (11) � � k EX � where the numerical wave number k ∆ is implicitly determined by the corresponding dispersion relation and k EX is the exact wave number for the given model. We wish to examine the phase error as a function of ω ∆ t in the range [0 , π ]. ∆ t is determined by h τ τ m , while ∆ x = ∆ y determined by CFL condition. We note that ω ∆ t = 2 π/ N ppp , where N ppp is the number of points per period, and is related to the number of points per wavelength as, N ppw = √ ǫ ∞ ν N ppp . We assume a uniform distribution and the following parameters which are appropriate constants for modeling aqueous Debye type materials: τ m = 8 . 1 × 10 − 12 sec , ǫ ∞ = 1 , ǫ s = 78 . 2 , τ r = 0 . 5 τ m . N. L. Gibson (Oregon State) Maxwell-PC Dispersive ICERM 2018 42 / 72
Stability and Dispersion Analyses PC-Debye FDTD Dispersion Analysis PC−Debye dispersion for FD with h τ =0.01, r=0.5 τ , θ =0 PC−Debye dispersion for FD with h τ =0.01, r=0.9 τ , θ =0 0.4 0.4 0.35 0.35 0.3 0.3 0.25 0.25 M=0 Φ 0.2 Φ 0.2 M=1 M=2 M=3 0.15 0.15 M=0 M=1 M=2 0.1 0.1 M=3 M=4 0.05 0.05 M=5 M=6 0 0 0 0.5 1 1.5 2 2.5 3 0 0.5 1 1.5 2 2.5 3 ω ∆ t ω ∆ t Figure: Plots of phase error at θ = 0 for (left column) τ r = 0 . 5 τ m , (right column) τ r = 0 . 9 τ m , using h τ = 0 . 01. N. L. Gibson (Oregon State) Maxwell-PC Dispersive ICERM 2018 43 / 72
Stability and Dispersion Analyses PC-Debye FDTD Dispersion Analysis PC−Debye dispersion for FD with h τ =0.001, r=0.5 τ , θ =0 PC−Debye dispersion for FD with h τ =0.001, r=0.9 τ , θ =0 0.5 0.5 0.45 0.45 0.4 0.4 0.35 0.35 0.3 0.3 M=0 M=1 Φ 0.25 Φ 0.25 M=2 M=3 0.2 0.2 M=4 M=5 0.15 0.15 M=6 M=0 0.1 0.1 M=1 M=2 0.05 0.05 M=3 0 0 0 0.5 1 1.5 2 2.5 3 0 0.5 1 1.5 2 2.5 3 ω ∆ t ω ∆ t Figure: Plots of phase error at θ = 0 for (left column) τ r = 0 . 5 τ m , (right column) τ r = 0 . 9 τ m , using h τ = 0 . 001. N. L. Gibson (Oregon State) Maxwell-PC Dispersive ICERM 2018 44 / 72
Stability and Dispersion Analyses PC-Debye FDTD Dispersion Analysis PC−Debye dispersion for FD with h τ =0.01, r=0.5 τ , ωτ µ =1 PC−Debye dispersion for FD with h τ =0.01, r=0.9 τ , ωτ µ =1 90 M=0 90 M=0 −1 −1 M=1 M=1 120 60 120 60 M=2 M=2 M=3 M=3 −3 −3 M=4 150 23 150 23 M=5 M=6 −5 −5 180 0 180 0 210 330 210 330 240 300 240 300 270 270 Figure: Log plots of phase error versus θ with fixed ω = 1 /τ m for (left column) τ r = 0 . 5 τ m , (right column) τ r = 0 . 9 τ m , using h τ = 0 . 01. Legend indicates degree M of the PC expansion. N. L. Gibson (Oregon State) Maxwell-PC Dispersive ICERM 2018 45 / 72
Stability and Dispersion Analyses PC-Debye FDTD Dispersion Analysis PC−Debye dispersion for FD with h τ =0.001, r=0.5 τ , ωτ µ =1 PC−Debye dispersion for FD with h τ =0.001, r=0.9 τ , ωτ µ =1 90 M=0 90 M=0 −1 −1 M=1 M=1 120 60 120 60 M=2 M=2 M=3 M=3 −3 −3 M=4 150 23 150 23 M=5 M=6 −5 −5 180 0 180 0 210 330 210 330 240 300 240 300 270 270 Figure: Log plots of phase error versus θ with fixed ω = 1 /τ m for (left column) τ r = 0 . 5 τ m , (right column) τ r = 0 . 9 τ m , using h τ = 0 . 001. Legend indicates degree M of the PC expansion. N. L. Gibson (Oregon State) Maxwell-PC Dispersive ICERM 2018 46 / 72
Random Debye Summary Summary We have presented a random ODE model for polydispersive Debye media N. L. Gibson (Oregon State) Maxwell-PC Dispersive ICERM 2018 47 / 72
Random Debye Summary Summary We have presented a random ODE model for polydispersive Debye media We described an efficient numerical method utilizing polynomial chaos (PC) and finite difference time domain (FDTD) N. L. Gibson (Oregon State) Maxwell-PC Dispersive ICERM 2018 47 / 72
Random Debye Summary Summary We have presented a random ODE model for polydispersive Debye media We described an efficient numerical method utilizing polynomial chaos (PC) and finite difference time domain (FDTD) Exponential convergence in the number of PC terms was demonstrated N. L. Gibson (Oregon State) Maxwell-PC Dispersive ICERM 2018 47 / 72
Random Debye Summary Summary We have presented a random ODE model for polydispersive Debye media We described an efficient numerical method utilizing polynomial chaos (PC) and finite difference time domain (FDTD) Exponential convergence in the number of PC terms was demonstrated We have proven (conditional) stability of the scheme via energy decay N. L. Gibson (Oregon State) Maxwell-PC Dispersive ICERM 2018 47 / 72
Random Debye Summary Summary We have presented a random ODE model for polydispersive Debye media We described an efficient numerical method utilizing polynomial chaos (PC) and finite difference time domain (FDTD) Exponential convergence in the number of PC terms was demonstrated We have proven (conditional) stability of the scheme via energy decay We have derived a discrete dispersion relation and computed phase errors N. L. Gibson (Oregon State) Maxwell-PC Dispersive ICERM 2018 47 / 72
Maxwell-Random Lorentz system Lorentz Model Lorentz Model We employ the physical assumption that electrons behave as damped harmonic oscillators, x + m ω 2 m ¨ x + 2 m ν ˙ 0 x = F driving . The polarization is then defined as electron dipole moment density: P + 2 ν ˙ ¨ P + ω 2 0 P = ǫ 0 ω 2 p E where ω 0 is the resonant frequency, ν is a damping coefficient, and ω p is referred to as a plasma frequency defined by ω 2 p = ( ǫ s − ǫ ∞ ) ω 2 0 . N. L. Gibson (Oregon State) Maxwell-PC Dispersive ICERM 2018 48 / 72
Maxwell-Random Lorentz system Lorentz Model Complex Permittivity Taking a Fourier transform of D = ǫ E + P and inserting the convolution form of the polarization model in for P , we get ˆ D ( ω ) = ǫ 0 ǫ ( ω ) ˆ E ( ω ) where ω 2 p ǫ ( ω ) = ǫ ∞ + 0 − ω 2 − i 2 νω . ω 2 For multiple Lorentz poles, the complex permittivity includes a (weighted) sum of mechanisms: N p ω 2 p , i � ǫ ( ω ) = ǫ ∞ + 0 , i − ω 2 − i 2 ν i ω. ω 2 i =1 N. L. Gibson (Oregon State) Maxwell-PC Dispersive ICERM 2018 49 / 72
Maxwell-Random Lorentz system Maxwell-Random Lorentz Random Polarization The multi-pole Lorentz model motivates a model with a continuum of Lorentz mechanisms, i.e., a distribution of dielectric parameters. We define a random polarization to be a function of a dielectric parameter treated as a random variable. The random Lorentz model is P + 2 ν ˙ ¨ P + ω 2 0 P = ǫ 0 ω 2 p E with parameter ω 2 0 treated as a random variable with probability distribution F on the interval ( a , b ). The macroscopic polarization is taken to be the expected value of the random polarization, � b P ( t , z ; ω 2 0 ) dF ( ω 2 P ( t , z ) = 0 ) . a N. L. Gibson (Oregon State) Maxwell-PC Dispersive ICERM 2018 50 / 72
Maxwell-Random Lorentz system Maxwell-Random Lorentz Random Polarization Figure: ω 2 0 ∼ U (0.75 ω 2 0 ,1.25 ω 2 0 ) N. L. Gibson (Oregon State) Maxwell-PC Dispersive ICERM 2018 51 / 72
Maxwell-Random Lorentz system Maxwell-Random Lorentz Complex Permittivity with random ω 2 0 Separate complex permittivity into real and imaginary parts ( ǫ = ǫ r + i ǫ i ): ω 2 p ( ω 2 0 − ω 2 ) ǫ r = ǫ ∞ + 0 − ω 2 ) 2 + 4 ν 2 ω 2 ( ω 2 2 ω 2 p νω ǫ i = 0 − ω 2 ) 2 + 4 ν 2 ω 2 . ( ω 2 Analytic integration is possible for uniform distribution: � b ω 2 1 b 0 ) 2 − 2 ω 2 0 ω 2 + ω 4 + 4 ν 2 ω 2 � � p ǫ r d ω 2 ln( ω 2 � E [ ǫ r ] = 0 = ǫ ∞ + � b − a 2( b − a ) � a a � b � ω 2 − ω 2 ω 2 1 � � b p ǫ i d ω 2 0 E [ ǫ i ] = 0 = ( b − a )arctan � b − a 2 νω � a a N. L. Gibson (Oregon State) Maxwell-PC Dispersive ICERM 2018 52 / 72
Maxwell-Random Lorentz system Maxwell-Random Lorentz Saltwater Data Figure: Fits for single-pole, saltwater data N. L. Gibson (Oregon State) Maxwell-PC Dispersive ICERM 2018 53 / 72
Maxwell-Random Lorentz system Maxwell-Random Lorentz Maxwell-Random Lorentz system In a polydisperse Lorentz material, we have ∂ E ∂ t = ∇ × H − ∂ P (14a) ǫ 0 ǫ ∞ ∂ t ∂ H ∂ t = − 1 ∇ × E (14b) µ 0 P + 2 ν ˙ ¨ P + ω 2 0 P = ǫ 0 ω 2 p E (14c) with � b P ( t , x ; ω 2 0 ) f ( ω 2 0 ) d ω 2 P ( t , x ) = 0 . a N. L. Gibson (Oregon State) Maxwell-PC Dispersive ICERM 2018 54 / 72
Maxwell-Random Lorentz system Maxwell-Random Lorentz Theorem (Stability of Maxwell-Random Lorentz) Let D ⊂ R 2 and suppose that E ∈ C (0 , T ; H 0 ( curl , D )) ∩ C 1 (0 , T ; ( L 2 ( D )) 2 ) , � 2 ) , and H ( t ) ∈ C 1 (0 , T ; L 2 ( D )) are solutions of P ∈ C 1 (0 , T ; L 2 (Ω) ⊗ L 2 ( D ) � the weak formulation for the Maxwell-Random Lorentz system along with PEC boundary conditions. Then the system exhibits energy decay E ( t ) ≤ E (0) ∀ t ≥ 0 , where the energy E ( t ) is defined as � √ ǫ 0 ǫ ∞ E ( t ) � √ µ 0 H ( t ) 2 2 ω 0 2 1 2 E ( t ) 2 = � � � � � � � � P ( t ) J ( t ) 2 + 2 + ω p √ ǫ 0 F + ω p √ ǫ 0 � � � � � � � � � � � � � � F (15) 2 ] and J := ∂ P where � u � 2 F = E [ � u � 2 ∂ t . Proof involves showing that � d E ( t ) = − 1 2 ν 2 � � J F ≤ 0 . � � dt E ( t ) ǫ 0 ω 2 � � p N. L. Gibson (Oregon State) Maxwell-PC Dispersive ICERM 2018 55 / 72
Maxwell-Random Lorentz system Maxwell-PC Lorentz Polynomial Chaos We wish to approximate the random polarization with orthogonal polynomials of the standard random variable ξ . Let ω 2 0 = r ξ + m and ξ ∈ [ − 1 , 1]. Suppressing the dimension of P and the spatial dependence, we have ∞ � α i ( t ) φ i ( ξ ) → ¨ P + 2 ν ˙ P + ω 2 0 P = ǫ 0 ω 2 P ( ξ, t ) = p E . i =0 Utilizing the Triple Recursion Relation for orthogonal polynomials: ξφ n ( ξ ) = a n φ n +1 ( ξ ) + b n φ n ( ξ ) + c n φ n − 1 ( ξ ) . the differential equation becomes ∞ � [ ¨ α i ( t ) + 2 ν ˙ α i ( t ) + m α i ( t )] φ i ( ξ ) i =0 ∞ � α i ( t ) [ a i φ i +1 ( ξ ) + b i φ i ( ξ ) + c i φ i − 1 ( ξ )] = ǫ 0 ω 2 + r p E φ 0 ( ξ ) . i =0 N. L. Gibson (Oregon State) Maxwell-PC Dispersive ICERM 2018 56 / 72
Maxwell-Random Lorentz system Maxwell-PC Lorentz Galerkin Projection We apply a Galerkin Projection onto the space of polynomials of degree at most p : α + 2 ν ˙ ¨ α = � � α + A � � f b 0 c 1 0 · · · 0 . . a 0 b 1 c 2 . ... ... ... A = rM + mI , M = . 0 0 . . . a p − 2 b b − 1 c p · · · 0 0 a p − 1 b p Or we can write as a first order system: α = � ˙ � β ˙ � α − 2 ν I � β + � β = − A � f , where � e 1 ǫ 0 ω 2 f = ˆ p E with ω p meaning expected value. N. L. Gibson (Oregon State) Maxwell-PC Dispersive ICERM 2018 57 / 72
Maxwell-Random Lorentz system Maxwell-PC Lorentz Maxwell-PC Lorentz The polynomial chaos system coupled with 1D Maxwell’s equations becomes ∂ E ∂ t = − ∂ H ǫ ∞ ǫ 0 ∂ z − β 0 ∂ H ∂ t = − 1 ∂ E µ 0 ∂ z α = � ˙ � β ˙ � α − 2 ν I � β + � β = − A � f Initial Conditions: α (0 , z ) = � E (0 , z ) = H (0 , z ) = � β (0 , z ) = 0 Boundary Conditions: E ( t , 0) = E L ( t ) and E ( t , z R ) = 0 N. L. Gibson (Oregon State) Maxwell-PC Dispersive ICERM 2018 58 / 72
Maxwell-Random Lorentz system Maxwell-PC Lorentz-FDTD We stagger three discrete meshes in the x and y directions and two discrete meshes in time: �� � � τ E x := | 0 ≤ ℓ ≤ L − 1 , 0 ≤ j ≤ J x ℓ + 1 2 , y j h �� � � τ E y | 0 ≤ ℓ ≤ L , 0 ≤ j ≤ J − 1 := x ℓ , y j + 1 h 2 �� � � τ H h := | 0 ≤ ℓ ≤ L − 1 , 0 ≤ j ≤ J − 1 x ℓ + 1 2 , y j + 1 2 τ E t := { ( t n ) | 0 ≤ n ≤ N } �� � � t n + 1 τ H t := | 0 ≤ n ≤ N − 1 . 2 N. L. Gibson (Oregon State) Maxwell-PC Dispersive ICERM 2018 59 / 72
Maxwell-Random Lorentz system Discrete Stability Staggered L 2 normed spaces Next, we define the L 2 normed spaces → R 2 | F = ( F x l + 1 � � h × τ E y F : τ E x 2 ) T , � F � E < ∞ − V E := 2 , j , F y l , j + 1 (18) h � � U : τ H V H := h − → R | U = ( U l + 1 2 ) , � U � H < ∞ (19) 2 , j + 1 with the following discrete norms and inner products L − 1 J − 1 � 2 , j | 2 + | F y ℓ, j + 1 2 | 2 � � � � F � 2 | F x ℓ + 1 , ∀ F ∈ V E E = ∆ x ∆ y (20) ℓ =0 j =0 L − 1 J − 1 � � � � ( F , G ) E = ∆ x ∆ y F x ℓ + 1 2 , j G x ℓ + 1 2 , j + F y ℓ, j + 1 2 G y ℓ, j + 1 , ∀ F , G ∈ V E (21) 2 ℓ =0 j =0 L − 1 J − 1 � � � U � 2 2 | 2 , ∀ U ∈ V H H = ∆ x ∆ y | U ℓ + 1 (22) 2 , j + 1 ℓ =0 j =0 L − 1 J − 1 � � ( U , V ) H = ∆ x ∆ y U ℓ + 1 2 V ℓ + 1 2 , ∀ U , V ∈ V H . (23) 2 , j + 1 2 , j + 1 ℓ =0 j =0 N. L. Gibson (Oregon State) Maxwell-PC Dispersive ICERM 2018 60 / 72
Maxwell-Random Lorentz system Discrete Stability We define a space and inner product for the random polarization in α and � β are now 2 × p + 1 matrices: vector notation, since � → R 2 × R p +1 � � h × τ E y α : τ E x − α = [ α 0 , . . . , α p ] , α k ∈ V E , � � V α := � � � α � h where the discrete L 2 grid norm and inner product are defined as p � α � 2 � α k � 2 � � α = E , ∀ � α ∈ V α k =0 p � � � α , � α , � ( � β ) α = ∀ � β ∈ V α . α k , β k E , k =0 We choose both spatial steps to be uniform and equal (∆ x = ∆ y = h ), and require that the usual CFL condition for two dimensions holds: √ 2 c ∞ ∆ t ≤ h . (24) N. L. Gibson (Oregon State) Maxwell-PC Dispersive ICERM 2018 61 / 72
Maxwell-Random Lorentz system Discrete Stability Theorem (Energy Decay for Maxwell-PC Lorentz-FDTD) If the stability condition (24) is satisfied, then the Yee scheme for the 2D TE mode Maxwell-PC Lorentz system satisfies the discrete identity 2 � � � n + 1 − 1 2 ν n + 1 � � 2 � δ t E = β (25) 2 � � h h n + 1 ǫ 0 ω 2 � � E 2 p � � A h for all n where 1 / 2 2 2 � � � � � � ω 2 2 ) H + �√ ǫ 0 ǫ ∞ E n � 2 1 � � � � µ 0 ( H n + 1 2 , H n − 1 � 0 E n α n β n h = E + � + � � � � ǫ 0 ω 2 ǫ 0 ω 2 � � � � p p � � � � α α (26) defines a discrete energy. α � 2 In the above � � A := ( A � α , � α ) α given A positive definite, which is true iff r < m . α � 2 α ≈ � E [ P ] � 2 2 + � StdDev( P ) � 2 2 = E [ �P� 2 2 ] = �P� 2 Note that � � F so that this is a natural extension of the Maxwell-Random Lorentz energy (15). N. L. Gibson (Oregon State) Maxwell-PC Dispersive ICERM 2018 62 / 72
Maxwell-Random Lorentz system Random Lorentz Dispersion Relation Theorem The dispersion relation for the Maxwell-Random Lorentz system is given by ω 2 c 2 ǫ ( ω ) = � k � 2 where the expected complex permittivity is given by ω 2 � � 0 ǫ ( ω ) = ǫ ∞ + ( ǫ s − ǫ ∞ ) E . 0 − ω 2 − i 2 νω ω 2 Where k is the wave vector and c = 1 / √ µ 0 ǫ 0 is the speed of light. The exact dispersion relation can be compared with a discrete dispersion relation to determine the amount of dispersion error. N. L. Gibson (Oregon State) Maxwell-PC Dispersive ICERM 2018 63 / 72
Maxwell-Random Lorentz system PC-Lorentz FDTD Dispersion Analysis Dispersion Error We define the phase error Φ for a scheme applied to a model to be � � k EX − k ∆ � � Φ = � , (27) � � k EX � where the numerical wave number k ∆ is implicitly determined by the corresponding discrete dispersion relation and k EX is the exact wave number for the given model. N. L. Gibson (Oregon State) Maxwell-PC Dispersive ICERM 2018 64 / 72
Maxwell-Random Lorentz system PC-Lorentz FDTD Dispersion Analysis Dispersion Error We define the phase error Φ for a scheme applied to a model to be � � k EX − k ∆ � � Φ = � , (27) � � k EX � where the numerical wave number k ∆ is implicitly determined by the corresponding discrete dispersion relation and k EX is the exact wave number for the given model. We wish to examine the phase error as a function of ω in the range around ω 0 . ∆ t is determined by h := ω 0 ∆ t / (2 π ), while ∆ x = ∆ y are determined by the CFL condition. We assume a uniform distribution and the following parameters Lorentz material: ν = 2 . 8 × 10 15 1/sec , ω 0 = 4 × 10 16 rad/sec . ǫ ∞ = 1 , ǫ s = 2 . 25 , N. L. Gibson (Oregon State) Maxwell-PC Dispersive ICERM 2018 64 / 72
Maxwell-Random Lorentz system PC-Lorentz FDTD Dispersion Analysis 0.7 PC Lorentz Dispersion with h=0.1 and p=1 2 r=0.1 0 2 r=0.2 0.6 0 2 r=0.3 0 2 0.5 r=0.4 0 2 r=0.5 0 0.4 0.3 0.2 0.1 0 0 2 4 6 8 10 16 Figure: Plots of phase error at θ = 0. N. L. Gibson (Oregon State) Maxwell-PC Dispersive ICERM 2018 65 / 72
Maxwell-Random Lorentz system PC-Lorentz FDTD Dispersion Analysis h=0.1 and r=0.1 w0^2 0.7 P=1 P=2 0.6 P=3 P=4 0.5 0.4 0.3 0.2 0.1 0 0 2 4 6 8 10 16 Figure: Plots of phase error at θ = 0. N. L. Gibson (Oregon State) Maxwell-PC Dispersive ICERM 2018 66 / 72
Maxwell-Random Lorentz system PC-Lorentz FDTD Dispersion Analysis h=0.1 and r=0.5 w0^2 0.7 P=1 P=2 0.6 P=3 P=4 0.5 0.4 0.3 0.2 0.1 0 0 2 4 6 8 10 16 Figure: Plots of phase error at θ = 0. N. L. Gibson (Oregon State) Maxwell-PC Dispersive ICERM 2018 67 / 72
Maxwell-Random Lorentz system PC-Lorentz FDTD Dispersion Analysis h=0.01 and r=0.1 w0^2 0.08 P=1 P=2 0.07 P=3 P=4 0.06 0.05 0.04 0.03 0.02 0.01 0 0 2 4 6 8 10 16 Figure: Plots of phase error at θ = 0. N. L. Gibson (Oregon State) Maxwell-PC Dispersive ICERM 2018 68 / 72
Maxwell-Random Lorentz system PC-Lorentz FDTD Dispersion Analysis h=0.01 and r=0.5 w0^2 0.7 P=1 P=2 0.6 P=3 P=4 P=5 0.5 P=6 0.4 0.3 0.2 0.1 0 0 2 4 6 8 10 16 Figure: Plots of phase error at θ = 0. N. L. Gibson (Oregon State) Maxwell-PC Dispersive ICERM 2018 69 / 72
Maxwell-Random Lorentz system PC-Lorentz FDTD Dispersion Analysis 10 -3 h=0.001 and r=0.1 w0^2 8 P=2 P=3 7 P=4 6 5 4 3 2 1 0 0 2 4 6 8 10 16 Figure: Plots of phase error at θ = 0. N. L. Gibson (Oregon State) Maxwell-PC Dispersive ICERM 2018 70 / 72
Future Directions Future Directions Extend to Drude meta-material models nonlinear polarization models viscoelastic system (partially done) Inverse problems from (actual) time-domain data N. L. Gibson (Oregon State) Maxwell-PC Dispersive ICERM 2018 71 / 72
References Bokil, V. A. & Gibson, N. L. (2012), Analysis of Spatial High Order Finite Difference Methods for Maxwell’s equations in dispersive media, IMA J. Numer. Anal. 32 (3): 926-956. Bokil, V. A. & Gibson, N. L. (2013), Stability and Dispersion Analysis of High Order FDTD Methods for Maxwell’s Equations in Dispersive Media, Contemporary Mathematics , Vol. 586. Bokil, V. A. & Gibson, N. L. (2014), Convergence Analysis of Yee Schemes for Maxwell’s Equations in Debye and Lorentz Dispersive Media, IJNAM , 11(4), 657-687. Gibson, N. L. (2015), A Polynomial Chaos Method for Dispersive Electromagnetics, Comm. in Comp. Phys. , vol. 18, issue 5, pp 1234-1263. Alvarez, Jacky & Fisher, Andrew (2017), Approximating Dispersive Materials With Parameter Distributions in the Lorentz Model, Oregon State University Math REU. N. L. Gibson (Oregon State) Maxwell-PC Dispersive ICERM 2018 72 / 72
Recommend
More recommend