Efficient method for evaluation of second-harmonic generation by surface integral equation

Lei Zhang; Shifei Tao; Zhenhong Fan; Rushan Chen

doi:10.1364/OE.25.028010

1. Introduction

Second-harmonic generation (SHG) is a nonlinear optical phenomenon constrained by material symmetry, where two photons of the same frequency interact with the material and generate one photon of twice the frequency [1]. In recent years, research of SHG from noble metal nanoparticles has drawn highly attention for its promising potential to probe physical and chemical properties of noble metal materials in the surface and interface. For non-centrosymmetric noble metal materials, SHG is considered to occur both in the bulk of the particle and on its surface [2,3]. However, surface contribution is sufficient to effectively represent the total SHG from gold structures, since it can be orders of magnitude greater than bulk contribution with both theoretical calculation and experimental results [4,5].

The nonlinear electromagnetic (EM) phenomena are usually modeled by time domain formulations [6,7]. However, it is difficult for time domain methods to deal with the dispersion problem without utilizing measured parameters of the material directly. Thus, frequency domain presents a rising trend to model and solve nonlinear optical processes. Jouni Mӓkitalo et al. [8] proposed a boundary element method (BEM) to evaluate the local-surface SHG from arbitrarily shaped particles. Then J. Butet et al. [9] extended BEM with the periodic Green’s function to calculate local-surface SHG from periodic structures. To take both nonlocal-bulk and local-surface contributions into consideration, Benedetti et al. [10,11] developed a less efficient formulation due to requirement of volume discretization based on volume integral equation (VIE). C. Forestiere et al. [12] calculated the SHG by only surface integral equation (SIE), taking both nonlocal-bulk and local-surface contributions into consideration with the equivalent surface electric and magnetic current. M. Luo et al. [13] developed a spectral element method (SEM) to solve the fundamental and second-harmonic (SH) fields simultaneously and self-consistently with only the bulk contribution. Xiaoyan Y. Z. Xiong et al. [14] took the mutual coupling between fundamental and SH field into consideration, solving the coupled-wave equations for second-harmonic generation by iteratively to design a compact nonlinear Yagi-Uda nanoantennas. It is worth mentioning that all the aforementioned methods are invoking the Love’s equivalent principle [15] during the SIE solving process.

In this express, we present a novel full wave analysis for evaluating surface SHG radiation with the tangential Poggio–Miller–Chang–Harrington–Wu–Tsai (PMCHWT) formulation in frequency domain. According to the new formulations, the EM fields at both fundamental frequency and SH, as well as the scattering cross section (SCS), are calculated by only half unknowns compared with the conventional SIE method, contributing to large conservation in computational memory and time. The proposed method also sharply declines the condition number of impendence matrix, and provides an efficient and promising approach for SHG evaluation with metal nanoparticles array and optimization design of nonlinear nanoantennas.

The remainder of the paper is organized as follows: Section 2 successively gives the detailed derivations of the conventional SIE method and our proposed SIE method for SHG evaluation after statement of the nonlinear problem. In section 3, some numerical results are presented to demonstrate the accuracy and efficiency of the proposed method. A brief conclusion with reference to future work is given in section 4.

2. Theory and formulations

Under the undepleted-pump approximation [16], the nonlinear problem of SHG evaluation in our implementation can be solved in three steps: 1) Solve a linear scattering problem with incident fields at the fundamental frequency. 2) Calculate the distribution of second-order nonlinear source polarization with the fundamental fields. 3) Solve a linear scattering problem at the SH frequency with specific SH incident fields obtained by the polarization induced by the fundamental harmonic field. This dealing is built on low SHG conversion efficiency, ignoring the mutual coupling between fundamental and SH field, which is justified by the fact that the measured SH signals are always orders of magnitude weaker than the source at the fundamental frequency [8].

2.1 The conventional SIE method

The solution domain of the EM field is divided into the interior domain $V_{i}$ , the exterior domain $V_{e}$ and the interface $S$ . The closed surface $S$ is assumed to be piece-wise smooth and oriented with its unitary normal $n$ points outward. We denote with $(E_{l}^{(ω)}, H_{l}^{(ω)})$ the fundamental total fields at angular frequency $ω$ , and with $(E_{l}^{(2 ω)}, H_{l}^{(2 ω)})$ the SH total fields at angular frequency $2 ω$ in the domain $V_{l}$ with $l = i, e$ . Furthermore, we denote with $(E_{0}^{(ω, i n c)}, H_{0}^{(ω, i n c)})$ the incident fields in $V_{e}$ , with $(ε_{i}^{(ν)}, μ_{i})$ the linear permittivity at frequency $ν = ω, 2 ω$ and the permeability of the interior domain, and with $(ε_{e}^{(Ω)}, μ_{e})$ the permittivity and the permeability of the exterior domain. Then we introduce the scattered fields $(E_{l}^{(ω, s c a)}, H_{l}^{(ω, s c a)})$ with $l = i, e$ at the fundamental frequency, defined as:

{\begin{matrix} E_{e}^{(ω, s c a)} = E_{e}^{(ω)} - E_{0}^{(ω, i n c)} \\ H_{e}^{(ω, s c a)} = H_{e}^{(ω)} - H_{0}^{(ω, i n c)} \end{matrix} i n V_{e} {\begin{matrix} E_{i}^{(ω, s c a)} = E_{i}^{(ω)} \\ H_{i}^{(ω, s c a)} = H_{i}^{(ω)} \end{matrix} i n V_{i}

The scattered fields at the fundamental frequency can be determined by solving the problem:

\begin{array}{l} {\begin{cases} \nabla \times E_{l}^{(ω, s c a)} = - j ω μ_{l} H_{l}^{(ω, s c a)} \\ \nabla \times H_{l}^{(ω, s c a)} = j ω ε_{l}^{(ω)} E_{l}^{(ω, s c a)} \end{cases} i n V_{l}, w i t h l = i, e, \\ {\begin{cases} n \times (E_{e}^{(ω, s c a)} - E_{i}^{(ω, s c a)}) = - n \times E_{0}^{(ω, i n c)} \\ n \times (H_{e}^{(ω, s c a)} - H_{i}^{(ω, s c a)}) = - n \times H_{0}^{(ω, i n c)} \end{cases} o n S \end{array}

with the radiation condition at infinity [17].

Once problem (2) has been solved, the local-surface SH polarization sources can be calculated. In this paper, only surface SHG is considered. The local-surface SH polarization field $P^{S}$ is decided by

P^{S} = ε_{0} {\overset{\leftrightarrow}{χ}}^{(2)} : E_{i}^{(ω)} E_{i}^{(ω)}

where

{\overset{\leftrightarrow}{χ}}^{(2)}

is the second-order surface nonlinear susceptibility,

E_{i}^{(ω)}

is the fundamental total field inside the interface S. The surface susceptibility tensor has only three non-vanishing and independent elements for centrosymmetric noble metal materials [1]. Equation (3) can be further written as

P^{S} = ε_{0} [χ_{⊥ ⊥ ⊥}^{(2)} n n n + χ_{⊥ | | | |}^{(2)} (n t_{1} t_{1} + n t_{2} t_{2}) + χ_{| | ⊥ | |}^{(2)} (t_{1} n t_{1} + t_{2} n t_{2})] : E_{i}^{(ω)} E_{i}^{(ω)}

where

(t_{1}, t_{2}, n)

is a system of three orthogonal vectors locally defined on the particle surface.

Once the local-surface second-harmonic polarization field $P^{S}$ has been obtained, the non-homogeneous problem of the SH fields is solved by

\begin{array}{l} {\begin{cases} \nabla \times E_{l}^{(2 ω)} = - j 2 ω μ_{l} H_{l}^{(2 ω)} \\ \nabla \times H_{l}^{(2 ω)} = j 2 ω ε_{l}^{(2 ω)} E_{l}^{(2 ω)} \end{cases} i n V_{l}, w i t h l = i, e, \\ {\begin{cases} n \times (E_{e}^{(2 ω)} - E_{i}^{(2 ω)}) = - j 2 ω P_{t}^{S} \\ n \times (H_{e}^{(2 ω)} - H_{i}^{(2 ω)}) = \frac{1}{ε^{'}} n \times \nabla_{S} P_{n}^{S} \end{cases} o n S \end{array}

with the radiation condition at infinity [17], where the surface electric and magnetic current are introduced with the contribution of tangential and normal component of the surface nonlinear polarization [18], i.e. $P_{t}^{S} = - n \times n \times P^{S},$ $P_{n}^{S} = P^{S} \cdot n$ , and

ε^{'}

is called the selvedge region permittivity [19].

In the conventional SIE method for SHG problem of Eqs. (2) and (5), the sources are removed by invoking Love’s equivalent principle [15]. The equivalent currents $(j_{e}^{(ν, e)}, j_{e}^{(ν, m)})$ with positioned outside the surface S, produce the field $(E_{e}^{(ν, s c a)}, H_{e}^{(ν, s c a)})$ in domain $V_{e}$ , and null-fields in domain $V_{i}$ , while the equivalent currents $(j_{i}^{(ν, e)}, j_{i}^{(ν, m)})$ defined inside the surface S, produce the original fields in domain $V_{i}$ and null fields in domain $V_{e}$ . By combining the tangential electric field integral equation (T-EFIE) and tangential magnetic field integral equation (T-MFIE) on both sides of the interface, the PMCHWT-2 formulation involved unknowns of two domains is obtained:

C_{2}^{(ν)} x_{2}^{(ν)} = y_{2}^{(ν)}

where the operator

C_{2}^{(ν)}

, the vector of unknowns

x_{2}^{(ν)}

, and the excitation vector

y_{2}^{(ν)}

are defined as

C_{2}^{(ν)} = (\begin{matrix} L_{e}^{(ν)} & K_{e}^{(ν)} & - L_{i}^{(ν)} & - K_{i}^{(ν)} \\ - K_{e}^{(ν)} & η_{e}^{- 2} L_{e}^{(ν)} & K_{i}^{(ν)} & - η_{i}^{- 2} L_{i}^{(ν)} \\ I & 0 & I & 0 \\ 0 & I & 0 & I \end{matrix}), x_{2}^{(ν)} = (\begin{matrix} j_{e}^{(ν, e)} \\ j_{e}^{(ν, m)} \\ j_{i}^{(ν, e)} \\ j_{i}^{(ν, m)} \end{matrix}), y_{2}^{(ν)} = (\begin{matrix} \frac{1}{2} n \times π_{S}^{(ν, m)} \\ - \frac{1}{2} n \times π_{S}^{(ν, e)} \\ π_{S}^{(ν, e)} \\ π_{S}^{(ν, m)} \end{matrix})

where the surface electric

π_{S}^{(ν, e)}

and magnetic

π_{S}^{(ν, m)}

current are defined as

{\begin{cases} π_{S}^{(ω, e)} = - n \times H_{0}^{(ω, i n c)} \\ π_{S}^{(ω, m)} = n \times E_{0}^{(ω, i n c)} \end{cases}, {\begin{cases} - π_{S}^{(2 ω, m)} = - j 2 ω P_{t}^{S} \\ π_{S}^{(2 ω, e)} = \frac{1}{ε^{'}} n \times \nabla_{S} P_{n}^{S} \end{cases}

with the definition of the

L_{l}^{(ν)}

and

K_{l}^{(ν)}

operators are given as follows

\begin{array}{l} L_{l}^{(ν)} {X (r)} = - n \times n \times i ν μ_{l} \int_{S} G_{l}^{(ν)} (r - r^{'}) X (r^{'}) d S^{'} \\ - n \times n \times \frac{1}{i ν ε_{l}} \int_{S} \nabla^{'} G_{l}^{(ν)} (r - r^{'}) {\nabla^{'}}_{S} \cdot X (r^{'}) d S^{'} \end{array}

K_{l}^{(ν)} {X (r)} = n \times n \times \int_{S} X (r^{'}) \times \nabla^{'} G_{l}^{(ν)} (r - r^{'}) d S^{'}

2.2. Our proposed SIE method

Different from the conventional SIE method introduced with the scattered fields, our proposed SIE method is built on the continuity of total fields across the interface. By this dealing, we obtain the PMCHWT-1 formulation for SHG evaluation at both fundamental frequency and second-harmonic involved unknowns of only one domain in interior or exterior. In this part, detailed derivations of difference between our proposed SIE method and the conventional SIE method will be shown. Simultaneously, we give the explanation of their consistency in theory.

The SIE applied by method of moments (MoM) to calculate the scattering field at the fundamental frequency according to Eq. (2) has been discussed in many literatures [20–23]. The total fields at the fundamental frequency can be determined by solving the problem:

\begin{array}{l} {\begin{cases} \nabla \times E_{l}^{(ω)} = - j ω μ_{l} H_{l}^{(ω)} \\ \nabla \times H_{l}^{(ω)} = j ω ε_{l}^{(ω)} E_{l}^{(ω)} \end{cases} i n V_{l}, w i t h l = i, e, \\ {\begin{cases} n \times (E_{e}^{(ω)} - E_{i}^{(ω)}) = n \times (E_{e}^{(ω, s c a)} + E_{0}^{(ω, i n c)} - E_{i}^{(ω, s c a)}) = 0 \\ n \times (H_{e}^{(ω)} - H_{i}^{(ω)}) = n \times (H_{e}^{(ω, s c a)} + H_{0}^{(ω, i n c)} - H_{i}^{(ω, s c a)}) = 0 \end{cases} o n S \end{array}

The boundary conditions of the total electric and magnetic fields on surface of nanoparticles at fundamental frequency can be indicated as:

n \times (E_{i}^{(ω)} - E_{e}^{(ω)}) = 0, n \times (H_{i}^{(ω)} - H_{e}^{(ω)}) = 0

where

n

is the outer normal of the interface. In our proposed method, the relationships of equivalent electric and magnetic currents

(J_{l}^{(ω)}, M_{l}^{(ω)})

and total fields

(E_{l}^{(ω)}, H_{l}^{(ω)})

in interior and external are given by:

{\begin{cases} J_{e}^{(ω)} = n_{1} \times H_{e}^{(ω)} \\ M_{e}^{(ω)} = - n_{1} \times E_{e}^{(ω)} \end{cases}, {\begin{cases} J_{i}^{(ω)} = n_{2} \times H_{i}^{(ω)} \\ M_{i}^{(ω)} = - n_{2} \times E_{i}^{(ω)} \end{cases}

Due to the only opposite direction of

n_{1}

and

n_{2}

, we can get the following with substituting Eq. (13) into Eq. (12):

J_{e}^{(ω)} = - J_{i}^{(ω)}, M_{e}^{(ω)} = - M_{i}^{(ω)}

Therefore, two groups of the unknown equivalent current $(J_{l}^{(ω)}, M_{l}^{(ω)})$ with $l = i, e$ at fundamental frequency can be reduced to only one group of interior or exterior. In our proposed SIE method, we use $(J^{(ω)}, M^{(ω)})$ to represent the equivalent electric and magnetic current at fundamental frequency for the sake of simplicity.

To compare similarity and difference between the conventional SIE method and our proposed SIE method, we take the solution of fundamental fields for example to explain. First, we list the equivalent electric and magnetic currents in the two SIE methods, as:

{\begin{matrix} j_{e}^{(ω, e)} = n \times H_{e}^{(ω, s c a)} \\ j_{e}^{(ω, m)} = - n \times E_{e}^{(ω, s c a)} \\ j_{i}^{(ω, e)} = - n \times H_{i}^{(ω, s c a)} \\ j_{i}^{(ω, m)} = n \times E_{i}^{(ω, s c a)} \end{matrix}, {\begin{matrix} J_{e}^{(ω)} = n \times H_{e}^{(ω)} \\ M_{e}^{(ω)} = - n \times E_{e}^{(ω)} \\ J_{i}^{(ω)} = - n \times H_{i}^{(ω)} \\ M_{i}^{(ω)} = n \times E_{i}^{(ω)} \end{matrix}

It can be observed that the equivalent electric and magnetic currents in the conventional SIE method are built on scattered fields while the equivalent electric and magnetic currents in our efficient SIE method are built on total fields. According to the relationship between scattered fields and total fields defined as Eq. (1), we can obtain:

{\begin{array}{l} J_{e}^{(ω)} = j_{e}^{(ω, e)} + n \times H_{0}^{(ω, i n c)} (a) \\ M_{e}^{(ω)} = j_{e}^{(ω, m)} - n \times E_{0}^{(ω, i n c)} (b) \\ J_{i}^{(ω)} = j_{i}^{(ω, e)} = - n \times H_{i}^{(ω, s c a)} (c) \\ M_{i}^{(ω)} = j_{i}^{(ω, m)} = n \times E_{i}^{(ω, s c a)} (d) \end{array}

As we can see from Eqs. (16.c) and (16.d), for two different SIE methods on SHG evaluation, the equivalent electric and magnetic currents in the interior of the metal domain are exactly the same. Thus, the same local-surface second-harmonic polarization field $P^{S}$ can be obtained by Eqs. (3) and (4), which is the key reason why our proposed SIE method achieve accurate SHG evaluation even if the mutual coupling between fundamental and second-harmonic field into consideration.

On the other hand, for two different SIE methods on SHG evaluation, the equivalent electric and magnetic currents in the exterior of the medium domain are not the same as we can see from Eqs. (16.a) and (16.b). However, different currents obtained by the two different SIE methods, $(j_{e}^{(ω, e)}, j_{e}^{(ω, m)})$ and $(J_{e}^{(ω)}, M_{e}^{(ω)})$ , will lead to the same scattered fields since the incident $(E_{0}^{(ω, i n c)}, H_{0}^{(ω, i n c)})$ have none contribution to the far-field. To validate this, we introduce $(J_{l}^{(ω, i n c)}, M_{l}^{(ω, i n c)})$ with $l = i, e$ to indicate the equivalent electric and magnetic currents in the interior and exterior domain. According to symmetry, we have:

J_{e}^{(ω, i n c)} = - J_{i}^{(ω, i n c)}, M_{e}^{(ω, i n c)} = - M_{i}^{(ω, i n c)}

Since the equivalent currents $(J_{i}^{(ω, i n c)}, M_{i}^{(ω, i n c)})$ produce the fields $(E_{0}^{(ω, i n c)}, H_{0}^{(ω, i n c)})$ in domain $V_{i}$ , and null scattered fields in domain $V_{e}$ , we can obtain:

{\begin{cases} E_{i \to e}^{(ω, s c a)} = L_{e}^{(ω)} (J_{i}^{(ω, i n c)}) + K_{e}^{(ω)} (M_{i}^{(ω, i n c)}) = 0 \\ H_{i \to e}^{(ω, s c a)} = - K_{e}^{(ω)} (J_{i}^{(ω, i n c)}) + η_{e}^{- 2} L_{e}^{(ω)} (M_{i}^{(ω, i n c)}) = 0 \end{cases}

where

(E_{i \to e}^{(ω, s c a)}, H_{i \to e}^{(ω, s c a)})

indicate scattered fields in domain

V_{e}

produced by

(J_{i}^{(ω, i n c)}, M_{i}^{(ω, i n c)})

. Due to the linearity of operators

(L_{l}^{(ν)}, K_{l}^{(ν)})

, we can also calculate scattered fields

(E_{e \to e}^{(ω, s c a)}, H_{e \to e}^{(ω, s c a)})

in domain

V_{e}

produced by

(J_{e}^{(ω, i n c)}, M_{e}^{(ω, i n c)})

, as:

{\begin{cases} E_{e \to e}^{(ω, s c a)} = L_{e}^{(ω)} (J_{e}^{(ω, i n c)}) + K_{e}^{(ω)} (M_{e}^{(ω, i n c)}) = 0 \\ H_{e \to e}^{(ω, s c a)} = - K_{e}^{(ω)} (J_{e}^{(ω, i n c)}) + η_{e}^{- 2} L_{e}^{(ω)} (M_{e}^{(ω, i n c)}) = 0 \end{cases}

Hence, the incident

(E_{0}^{(ω, i n c)}, H_{0}^{(ω, i n c)})

have none contribution to the far-field, which is the key reason why our efficient SIE method reduces the number of unknowns to only one group equivalent currents of interior or exterior.

In a word, the two SIE methods differ in the form, but keep the same in essence. Therefore, according to the similar solutions of Eqs. (2) and (5), the total fields at the SH can be analogically determined by solving the problem:

\begin{array}{l} {\begin{cases} \nabla \times E_{l}^{(2 ω)} = - j ω μ_{l} H_{l}^{(2 ω)} \\ \nabla \times H_{l}^{(2 ω)} = j ω ε_{l}^{(2 ω)} E_{l}^{(2 ω)} \end{cases} i n V_{l}, w i t h l = i, e, \\ {\begin{cases} n \times (E_{e}^{(2 ω)} - E_{i}^{(2 ω)}) = n \times (E_{e}^{(2 ω, s c a)} + E_{r}^{(2 ω, i n c)} - E_{i}^{(2 ω, s c a)}) = 0 \\ n \times (H_{e}^{(2 ω)} - H_{i}^{(2 ω)}) = n \times (H_{e}^{(2 ω, s c a)} + H_{r}^{(2 ω, i n c)} - H_{i}^{(2 ω, s c a)}) = 0 \end{cases} o n S \end{array}

where

(E_{r}^{(2 ω, i n c)}, H_{r}^{(2 ω, i n c)})

are the equivalent SH incident electric and magnetic fields obtained by the SH surface polarization source based on Maxwell’s equations as:

{\begin{cases} E_{r}^{(2 ω, i n c)} = (L_{i}^{(ω)} + L_{e}^{(ω)}) π_{S}^{(2 ω, e)} + (K_{i}^{(ω)} + K_{e}^{(ω)}) π_{S}^{(2 ω, m)} \\ H_{r}^{(2 ω, i n c)} = (- K_{i}^{(ω)} - K_{e}^{(ω)}) π_{S}^{(2 ω, e)} + (η_{i}^{- 2} L_{i}^{(ω)} + η_{e}^{- 2} L_{e}^{(ω)}) π_{S}^{(2 ω, m)} \end{cases}

where the subscript r emphasizes the equivalent incident wave is a non-uniform wave, differs from that of the fundamental frequency.

By combining the tangential electric field integral equation (T-EFIE) and tangential magnetic field integral equation (T-MFIE) on both sides of the interface, the PMCHWT-1 formulation involved unknowns of only one domain in exterior is obtained:

C_{1}^{(ν)} x_{1}^{(ν)} = y_{1}^{(ν)}

where the operator

C_{1}^{(ν)}

and the vector of unknowns

x_{1}^{(ν)}

are defined as

C_{1}^{(ν)} = (\begin{matrix} L_{e}^{(ν)} + L_{i}^{(ν)} & K_{e}^{(ν)} + K_{i}^{(ν)} \\ - K_{e}^{(ν)} - K_{i}^{(ν)} & η_{e}^{- 2} L_{e}^{(ν)} + η_{i}^{- 2} L_{i}^{(ν)} \end{matrix}), x_{1}^{(ν)} = (\begin{matrix} J^{(ν)} \\ M^{(ν)} \end{matrix}),

and the excitation vector

y^{(ν)}

is defined as

y_{1}^{(ω)} = (\begin{matrix} n \times n \times E_{0}^{(ω, i n c)} \\ n \times n \times H_{0}^{(ω, i n c)} \end{matrix}), y_{1}^{(2 ω)} = (\begin{matrix} n \times n \times E_{r}^{(2 ω, i n c)} \\ n \times n \times H_{r}^{(2 ω, i n c)} \end{matrix})

Hence, our SIE method for SHG evaluation not only reduces the number of unknowns, but also greatly shifts the impendence matrix blocks compared to the conventional SIE method, leading to remarkable computational effect.

In our implementation, all the integrals in the presented formulations can be calculated with high precision by utilizing the singularity subtraction technique [24,25] for the singular part of the inner products and high-order Gaussian quadrature [26] for the non-singular part of the inner products with smooth kernels.

2.4 Calculation of scattering cross sections

To investigate the SHG of nanoparticles illuminated by strong laser, the physical quantity of scattering cross sections (SCS) are defined as the SH radiated power with

σ = \lim_{r \to \infty} \frac{1}{4 π r^{2}} \frac{| E^{s} |^{2}}{| E^{i} |^{2}}

where

E^{s}

is the scattering field determined by the distribution of induced currents, and

E^{i}

is the incident field decided by laser source, respectively. The SCS is a significant characterization of energy transmission from near incident field to far field.

3. Numerical results

After formulations given, we start to model SHG from noble metal nanoparticles with arbitrary shapes. First, we give an example of spherical of radius r = 50nm to validate the accuracy and efficiency of our proposed SIE method. Then spherical particles of different radius are modeled to demonstrate applicability without a lack of generality. At last, we apply the method into SHG from a studied experimentally particle of L-shaped. The material of all models is gold and the exterior domain is assumed as vacuum. This yields $ε^{'} = ε_{0}$ . The related index of gold at optical frequency is obtained from [27].

3.1 The spherical particle of radius r = 50 nm

Consider a gold sphere with radius r = 50 nm excited by a uniform plane wave of unitary intensity, which is propagated along the $+ z$ axis and linearly polarized along x. Corresponding to the resonant frequency of the gold spherical particle, the incident wavelength is set as λ = 520 nm. The gold spherical particle is discretized with RWG basis [28] of N = 930.

Figure 1(a) plots the E-plane ( $φ = 0^{\circ}$ ) and H-plane ( $φ = 90^{\circ}$ ) SCS at the fundamental frequency. It can be observed the results calculated by the proposed SIE method and the Mie solutions [29] keep good agreement. Figures 1(b) and 1(c) show the normalized SCS at the SH frequency. In Fig. 1(b), only the normal components of local-surface polarization $P_{n}^{S}$ are active, while in Fig. 1(c) only the tangential components of local-surface polarization $P_{t}^{S}$ are considered. It can be found very good agreement between the proposed SH-SIE method and the SH-Mie formulations in all cases. Comparison of the convergence property for the iterative solver between the proposed SH-SIE method and the conventional SIE method is shown in Fig. 1(d). It is observed our method presents a remarkable convergence behavior with only half unknowns compared to the conventional SIE method.

Fig. 1 The results of spherical gold particle of radius r = 50nm. (a) SCS at fundament frequency. SCS at SH frequency with (b) $(χ_{⊥ ⊥ ⊥}, χ_{| | ⊥ | |}) = (1, 0)$ , (c) $(χ_{⊥ ⊥ ⊥}, χ_{| | ⊥ | |}) = (0, 1)$ . (d) Comparison of convergence property.

Download Full Size | PDF

Table 1 shows the comparison of CPU computational time and memory between the proposed SH-SIE method and the conventional SIE method. For specified model in same mesh, the proposed SH-SIE method not only presents high efficiency on CPU computational time, but also has capability to save much CPU computational memory.

Table 1. Comparison of Two Methods for SHG of Spherical Nanoparticles r = 50nm

View Table | View all tables in this article

3.2 The spherical particles of different radius

We next model spherical particles of different radius to demonstrate applicability without a lack of generality. The radius of spherical particles is chosen to be 10nm, 20nm and 100nm, respectively. Different methods are applied and different components of nonlinear source polarization are considered in SHG calculation. The results are shown in Fig. 2.

Fig. 2 The results of spherical gold particle with radius (a)-(c) 10nm, (d)-(f) 20nm, (g)-(i) 100nm. SCS at SH frequency with (b)(d)(g) $(χ_{⊥ ⊥ ⊥}, χ_{| | ⊥ | |}) = (1, 0)$ , (c)(e)(h) $(χ_{⊥ ⊥ ⊥}, χ_{| | ⊥ | |}) = (0, 1)$ . (d)(f)(i) Comparison of convergence property.

Download Full Size | PDF

In Figs. 2(a), 2(d) and 2(g), only the normal components of local-surface polarization $P_{n}^{S}$ are active, whereas in Figs. 2(b), 2(e) and 2(h) only the tangential components of local-surface polarization $P_{t}^{S}$ are considered for different radius r = 10nm, 20nm and 100nm, successively. It is observed that the results of proposed method are surprisingly accurate compared with SH-Mie formulations in all test cases. Figures 2(c), 2(f) and 2(i) show the convergence behavior of the iterative solver comparison. Similar remarkable convergence behavior is observed during comparison between the proposed SH-SIE method and the conventional SIE method.

Table 2 shows the comparison of CPU computational time and memory between the proposed SH-SIE method and the conventional SIE method of the three examples above. In all test cases, the proposed method presents high efficiency on CPU computational time and save much CPU computational memory. Furthermore, the more unknowns are related, the more prominent effects are performed.

Table 2. Comparison of Two Methods for SHG of Spherical Nanoparticles r = 10nm, 20nm and 100nm

View Table | View all tables in this article

3.3 The L-shaped particle

To further research the distribution of surface nonlinear polarization and properties of the SH radiation, we apply our method to a noncentrosymmetric L-shaped gold nanoparticle, which has been studied experimentally [30–33]. Figure 3 shows the gold L-shaped particle with the coordinate system. The particle consists of two axisymmetric arms of 150 nm × 50 nm with height of 20 nm. The single gold L-shaped nanoparticle is positioned in vacuum excited by a uniform plane wave of unitary intensity propagating in z-direction at its corresponding plasmon resonances. The effective second-order susceptibility components are $(χ_{⊥ ⊥ ⊥}, χ_{| | ⊥ | |}) = (250, 1)$ [5].

Fig. 3 The rounded-edge L-shaped particle.

Download Full Size | PDF

The amplitude of nonlinear surface polarization from the L-shaped gold nanoparticle by the incident wave at λx = 995 nm for x-polarization is shown as Fig. 4(a), and λy = 662 nm for y-polarization is shown as Fig. 4(b). These resonant plasmon wavelengths are decided by sensitive association with the charge oscillations along the arms and the arm length. It is obvious to see that the nonlinear source polarization most occurred at the corners of the particle, which is agreed with theoretical calculation by BEM in Ref [8]. and experimental measurements in supporting information in Ref [31]. Compared to polarization formation of the whole particle, only a small portion of the particle surface plays a significant role in giving rise to SHG. In another word, even a tiny flaw in the corner is likely to bring a sharp change in SHG. That’s why we model a rounded-edge L-shaped particle and have fine discretization at the corners to ensure sufficiently smooth around sharp edges and corners.

Fig. 4 The nonlinear surface polarization (a)(b) and SH radiation pattern (c)(d) from an L-shaped gold nanoparticles, with (a)(c) corresponds to incident wave with x-polarization and (b)(d) corresponds to incident wave with y-polarization.

Download Full Size | PDF

The full radiation patterns computed with the proposed method from the L-shaped gold nanoparticle by the incident wave at λx = 995 nm for x-polarization is shown as Fig. 4(c), and λy = 662 nm for y-polarization is shown as Fig. 4(d). To study the scattering field of second harmonics, which is the main concern in the L-shaped particle research, we should pay attention to the localization effect of nonlinear polarization source. Owing to the irregularity of this L-shaped nanoparticle, the local field is strongly enhanced with highly inhomogeneous distribution in local space, leading to in corresponding hot spots. Thus, it has to be considered that the signal arising from the polarization source could be strongly enhanced for the hot spots in the model of SH generation [5,34].When the incident plane wave is x-polarized, two hot spots arise at the end corner of arms of the L-shaped nanoparticle as shown in Fig. 4(a), whereas, in correspondence of a y-polarized excitation, another hot spot appears in the remaining corner of the L-shaped nanoparticle as shown in Fig. 4(b). As we can see, the most evident distinction caused by the new additional local field enhancement or hot spot is symmetric or non-symmetric between Figs. 4(c) and 4(d). Furthermore, according to Eq. (3), the shape of SH radiation pattern is affected by the distribution of the electric field on the surface of the particle at the fundamental frequency. For qualitative comparison, the highest power per goes to forward and backward directions for the case of x-polarized input, but not y-polarized input, which is also agreed with theoretical calculation by BEM in Ref [8]. To quantitatively verify our proposed method, we take the symmetry of highest power per unit to forward and backward directions into consideration. The results show that the relative errors are respectively 0.68% and 0.89% for the case of Figs. 4(c) and 4(d).

To find the experimental measurement data supporting our explanation, we compare the shape and feature of the scattering pattern with color graphs in the previous literature. In the experimental measurements, one is usually interested in the x- and y-components of the SH radiation in the forward direction. Symmetric considerations require that, in the case of an ideal particle, only the y-component can be non-zero [8,31,32], which is also observed in our SH radiation pattern. For validation purposes, we ensure that our results of SHG radiation by the proposed method following this symmetric rule.

On the premise of accuracy,our proposed SH-SIE method only need half unknowns compared with the conventional SIE method for the same model in same mesh. In theory, it is easy to understand half unknowns mean a quarter of computational memory. In addition, both the calculation of extra two terms at right hand and zero elements in impendence matrix lead to memory consumption [8,14]. On the other hand, the computational complexity of SIE method is of the order O(N³). It is predictable the proposed method is more efficient with the increasing of unknowns. More importantly, the new combination of impendence matrix blocks by the proposed method totally improves the condition number essentially.

4. Conclusion

A surface integral equation method is present as a full wave analysis to evaluate surface SHG from noble metal nanoparticles of virtually arbitrary shape. By introducing specific SH incident fields with the local SH surface polarization, only half unknowns of the equivalent electric and magnetic currents in the interior or exterior built on total fields are involved compared with the conventional SIE method.

Numerical examples of nanoparticles in different radius are given to demonstrate accuracy and efficiency. In all test cases, the proposed SIE method presents high efficiency on CPU computational time and save much CPU computational memory.

To further research the distribution of surface nonlinear polarization and properties of the SH radiation, we apply our method to a noncentrosymmetric L-shaped gold nanoparticle, which has been studied experimentally. The calculations suggest the nonlinear source polarization most occurred at the corners of the particle, which is justified by the experimental fact. Only a small portion of the particle surface plays a significant role in giving rise to SHG, so even a tiny flaw in the corner is likely to bring a sharp change in SHG.

At last, the proposed SIE method enables accurate and efficient simulation for nonlinear optics. This provides an efficient and promising approach for evaluation of nonlinear optical radiation generated from metal nanoparticles array and design of plasmonic metamaterials with special nonlinear optical properties in nano-scale.

Funding

This work was supported in part by Natural Science Foundation of 61431006, 61371037 and Defense Industrial Technology Development Program.

Acknowledgments

We thank Xiaoyan Y.Z. Xiong from the Department of Electrical and Electronic Engineering, the University of Hong Kong for helpful discussions about theoretical details and Wei E.I. Sha from Zhenjiang University for significant analysis of theoretical calculation.

References and links

1. Y. R. Shen, The principle of Nonlinear Optics (John Wiley & Sons Inc.,1984).

2. M. Kauranen and A. V. Zayats, “Nonlinear plasmonics,” Nat. Photonics 6, 737–748 (2012).

3. J. E. Sipe, V. Mizrahi, and G. I. Stegeman, “Fundamental difficulty in the use of second-harmonic generation as a strictly surface probe,” Phys. Rev. B Condens. Matter 35(17), 9091–9094 (1987). [PubMed]

4. P. Guyot-Sionnest, W. Chen, and Y. R. Shen, “General considerations on optical second-harmonic generation from surfaces and interfaces,” Phys. Rev. B Condens. Matter 33(12), 8254–8263 (1986). [PubMed]

5. F. Wang, F. Rodríguez, W. Albers, R. Ahorinta, J. Sipe, and M. Kauranen, “Surface and bulk contributions to the second-order nonlinear optical response of a gold film,” Phys. Rev. B 80(23), 308–310 (2009).

6. C. M. Dissanayake, M. Premaratne, I. D. Rukhlenko, and G. P. Agrawal, “FDTD modeling of anisotropic nonlinear optical phenomena in silicon waveguides,” Opt. Express 18(20), 21427–21448 (2010). [PubMed]

7. M. A. Alsunaidi, H. M. Al-Mudhaffar, and H. M. Masoudi, “Vectorial FDTD technique for the analysis of optical second-harmonic generation,” IEEE Photonics Technol. Lett. 21, 310–312 (2009).

8. J. Mäkitalo, S. Suuriniemi, and M. Kauranen, “Boundary element method for surface nonlinear optics of nanoparticles,” Opt. Express 19(23), 23386–23399 (2011). [PubMed]

9. J. Butet, B. Gallinet, K. Thyagarajan, and O. J. F. Martin, “Second harmonic generation from periodic arrays of arbitrary shape plasmonic nanostructures: a surface integral approach,” J. Opt. Soc. Am. B 30, 2970–2979 (2013).

10. A. Benedetti, M. Centini, C. Sibilia, and M. Bertolotti, “Engineering the second harmonic generation pattern from coupled gold nanowires,” J. Opt. Soc. Am. B 27, 408–416 (2010).

11. A. Benedetti, M. Centini, M. Bertolotti, and C. Sibilia, “Second harmonic generation from 3D nanoantennas: on the surface and bulk contributions by far-field pattern analysis,” Opt. Express 19(27), 26752–26767 (2011). [PubMed]

12. C. Forestiere, A. Capretti, and G. Miano, “Surface integral method for second harmonic generation in metal nanoparticles including both local-surface and nonlocal-bulk sources,” J. Opt. Soc. Am. B 30(9), 2355–2364 (2013).

13. M. Luo and Q. H. Liu, “Enhancement of second-harmonic generation in an air-bridge photonic crystal slab: simulation by spectral element method,” J. Opt. Soc. Am. B 28, 2879–2887 (2011).

14. X. Y. Z. Xiong, L. J. Jiang, W. E. I. Sha, Y. H. Lo, and W. C. Chew, “Compact Nonlinear Yagi-Uda Nanoantennas,” Sci. Rep. 6, 18872 (2016). [PubMed]

15. J. Jackson, Classical Electrodynamics (John Wiley & Sons, Inc., 1999).

16. R. W. Boyd, Nonlinear Optics, Third Edition (Academic, 2008).

17. J. Jin, Theory and Computation of Electromagnetic Fields (Wiley, 2010).

18. T. F. Heinz, Nonlinear Surface Electromagnetic Phenomena, H.-E. Ponath and G. I. Stegeman, eds. (Elsevier, 1991) Chapter 5, p. 397–405.

19. J. E. Sipe, V. C. Y. So, M. Fukui, and G. I. Stegeman, “Analysis of second-harmonic generation at metal surfaces,” Phys. Rev. B 21, 4389–4403 (1980).

20. P. Ylä-Oijala, M. Taskinen, and J. Sarvas, “Surface Integral Equation Method for General Composite Metallic and Dielectric Structures with Junctions,” Prog. Electromagnetics Res. 52(3), 81–108 (2005).

21. P. Ylä-Oijala, M. Taskinen, and S. Järvenpää, “Surface integral equation formulations for solving electromagnetic scattering problems with iterative methods,” Radio Sci. 40(6), 5391–5392 (2005).

22. Z. G. Qian and C. C. Weng, “Fast Full-Wave Surface Integral Equation Solver for Multiscale Structure Modeling,” IEEE Trans. Antenn. Propag. 57(11), 3594–3601 (2009).

23. P. Yla-Oijala and M. Taskinen, “Application of combined field Integral equation for electromagnetic scattering by dielectric and composite objects,” IEEE Trans. Antenn. Propag. 53(3), 1168–1173 (2005).

24. R. D. Graglia, “On the numerical integration of the linear shape functions times the 3-D Green’s function or its gradient on a plane triangle,” IEEE Trans. Antenn. Propag. 41(10), 1448–1455 (1993).

25. P. Yla-Oijala and M. Taskinen, “Calculation of CFIE impedance matrix elements with RWG and n×RWG functions,” IEEE Trans. Antenn. Propag. 51(8), 1837–1846 (2003).

26. D. A. Dunavant, “High degree efficient symmetrical Gaussian quadrature rules for the triangle,” Int. J. Numer. Methods Eng. 21(6), 1129–1148 (1985).

27. P. Johnson and R. Christy, “Optical constants of the noble metals,” Phys. Rev. B 6, 4370–4379 (1972).

28. M. S. Rao, D. R. Wilton, and A. W. Glisson, “Electromagnetic scattering by surfaces of arbitrary shape,” IEEE Trans. Antenn. Propag. 30(3), 409–418 (1982).

29. Y. Pavlyukh and W. Hubner, “Nonlinear Mie scattering from spherical particles,” Phys. Rev. B 70, 245–434 (2004).

30. S. Kujala, B. K. Canfield, M. Kauranen, Y. Svirko, and J. Turunen, “Multipole interference in the second-harmonic optical radiation from gold nanoparticles,” Phys. Rev. Lett. 98(16), 167403 (2007). [PubMed]

31. H. Husu, R. Siikanen, J. Mäkitalo, J. Lehtolahti, J. Laukkanen, M. Kuittinen, and M. Kauranen, “Metamaterials with tailored nonlinear optical response,” Nano Lett. 12(2), 673–677 (2012). [PubMed]

32. B. Canfield, S. Kujala, K. Jefimovs, J. Turunen, and M. Kauranen, “Linear and nonlinear optical responses influenced by broken symmetry in an array of gold nanoparticles,” Opt. Express 12(22), 5418–5423 (2004). [PubMed]

33. B. Lambrecht, A. Leitner, and F. R. Aussenegg, “Femtosecond decay-time measurement of electron-plasma oscillation in nanolithographically designed silver particles,” Appl. Phys. B 64(2), 269–272 (1997).

34. A. Benedetti, M. Centini, M. Bertolotti, and C. Sibilia, “Second harmonic generation from 3D nanoantennas: on the surface and bulk contributions by far-field pattern analysis,” Opt. Express 19(27), 26752–26767 (2011). [PubMed]

Method	Unknowns	Time(s)	Memory(MB)
The proposed SH-SIE method	7,704	389.05	46.86
The conventional SIE method	15,408	2889.92	287.36

Model	The proposed SH-SIE method			The conventional SIE method
Model	Unknowns	Time(s)	Memory(MB)	Unknowns	Time(s)	Memory(MB)
r = 10nm	1,860	61.63	5.260	3720	230.92	22.176
r = 20nm	4,614	178.79	27.128	9228	1139.31	696.140
r = 100nm	16,146	2452.03	1,060.632	32292	17688.42	8,178.304

Method	Unknowns	Time(s)	Memory(MB)
The proposed SH-SIE method	7,704	389.05	46.86
The conventional SIE method	15,408	2889.92	287.36

Model	The proposed SH-SIE method			The conventional SIE method
Model	Unknowns	Time(s)	Memory(MB)	Unknowns	Time(s)	Memory(MB)
r = 10nm	1,860	61.63	5.260	3720	230.92	22.176
r = 20nm	4,614	178.79	27.128	9228	1139.31	696.140
r = 100nm	16,146	2452.03	1,060.632	32292	17688.42	8,178.304

Efficient method for evaluation of second-harmonic generation by surface integral equation

Abstract

1. Introduction

2. Theory and formulations

2.1 The conventional SIE method

2.2. Our proposed SIE method

2.4 Calculation of scattering cross sections

3. Numerical results

3.1 The spherical particle of radius r = 50 nm

3.2 The spherical particles of different radius

3.3 The L-shaped particle

4. Conclusion

Funding

Acknowledgments

References and links

Cited By

Figures (4)

Tables (2)

Equations (25)

Optics Express