Accurate en efficiënte simulatie van nanoschaal optische vanging met FETI-DP en Maxwell-spanningstensor

Samenvatting

Numerieke simulatie speelt een belangrijke rol voor de voorspelling van optische vanging gebaseerd op plasmonische nano-optische pincetten. Echter, gecompliceerde structuren en drastische lokale veldversterking van plasmonische effecten brengen grote uitdagingen voor traditionele numerieke methoden. In dit artikel wordt een accurate en efficiënte numerieke simulatiemethode voorgesteld gebaseerd op een dual-primal finite element tearing and interconnecting (FETI-DP) en Maxwell-spanningstensor, om de optische kracht en potentieel voor het vangen van nanodeeltjes te berekenen. Een lage-rang sparsificatie benadering wordt geïntroduceerd om de FETI-DP simulatieprestaties verder te verbeteren. De voorgestelde methode kan een grootschalig en complex probleem ontbinden in kleinschalige en eenvoudige problemen door gebruik te maken van non-overlappende domeindeeling en flexibele mesh-discretisatie, wat hoge efficiëntie en paralleliseerbaarheid vertoont. Numerieke resultaten tonen de effectiviteit van de voorgestelde methode voor de voorspelling en analyse van optische vanging op nanoschaal.

Inleiding

Plasmonische optische pincetten gebaseerd op oppervlakteplasmonen (SPs) trekken veel aandacht en zijn breed toegepast om nanodeeltjes te vangen [1-6]. SP is een resonantieverschijnsel veroorzaakt door de koppeling van invallicht met een specifieke golfengte en vrije elektronen aan de interface van metalen en dielectrica [7]. SPs stellen de optische pincetten in staat om de diffractielimiet te doorbreken. Bovendien kan de drastische lokale veldversterking van de SPs de eis aan intensiteit van invallicht verminderen [7, 8]. Echter, SPs zijn nauw gerelateerd aan het materiaal en de afmetingen van objecten alsook de golfengte van invallicht, wat een groot aantal experimenten vereist om de optimale parameters van SP optische pincetten in de praktijk te bepalen. Op basis hiervan speelt de simulatiemethode een steeds belangrijker rol als hulpmiddel voor het ontwerp en de optimalisatie van SP optische pincetten [9]. In deze simulaties is de berekening van optische kracht vereist om een stabiele vanging te voorspellen. Voor regelmatige objecten zoals spheren kan de optische kracht analytisch worden afgeleid uit de gegeneraliseerde Lorenz-Mie theorie [10, 11]. Echter, voor objecten met gecompliceerde configuraties zijn numerieke methoden die de goberende Maxwellvergelijkingen rigoros oplossen noodzakelijk voor het modelleren van de elektromagnetische velden en de daaropvolgende optische kracht en potentieel.

Deze numerieke methoden kunnen hoofdzakelijk worden gecategoriseerd in differentiaalvergelijkingsmethoden (DEMs) en integraalvergelijkingsmethoden (IEMs) [12-15]. Vergeleken met IEMs tonen differentiaalvergelijkingsmethoden (DEMs) superieure capaciteiten in het omgaan met gecompliceerde geometrieën en componenten. DEMs hebben ook het voordeel van een directe berekening van nabijvelddistributie, wat een belangrijke rol speelt in de analyse van SPs. Als representatieve DEM is de finite-difference time-domain (FDTD) methode geïmplementeerd in het tijdsdomein, wat makkelijk breedbandinformatie en transeerresponsen kan verkrijgen [16, 17]. Echter, de FDTD eist een accurate dispersiemodel om de frequentie-afhankelijke materiaaleigenschappen in SPs te beschrijven, terwijl de FDTD oplossingsnauwkeurigheid sterk afhangt van de benaderingsnauwkeurigheid van dit dispersiemodel [18]. Daarnaast steunt de FDTD op gestructureerde meshes, wat vaak leidt tot trap-fouten voor gekromde oppervlakken. Als een andere representatieve DEM is de finite element method (FEM) breed geadopteerd omdat het gemakkelijk dispersief materiaal in het frequentiedomein kan hanteren en de trap-fout kan elimineren door ongestructureerde mesh [19-22]. Vergeleken met de FDTD kan de FEM direct gemeten materiaalparameters aan zonder enig benaderend dispersiemodel. Echter, drastische lokale veldversterkingen in de SPs vereisen fijnere meshes in de FEM-discretisatie. Bovendien zullen objecten met grote afmetingen en meerdere objecten het aantal onbekenden dramatisch verhogen. Deze factoren zullen slecht-geconditioneerde matricesystemen en enorme computationele consumpties veroorzaken, wat grote uitdagingen brengt voor traditionele FEM voor de analyse van SP-versterkte optische vanging.

In dit artikel wordt een efficiënte dual-primal finite element tearing and interconnecting (FETI-DP) methode geïntroduceerd om de optische vanging op nanoschaal te simuleren. De FETI-DP adopteert een non-overlappende domeindecompositie schema, wat een oorspronkelijk grootschalig complex probleem verdeelt in een reeks kleinschalige eenvoudige problemen om ze te overwinnen. Het dwingt een transmissievoorwaarde op de subdomeininterfaces af om de veldcontinuïteit te garanderen, en introduceert een duale variabele om het oorspronkelijke driedimensionale (3D) probleem te reduceren tot een tweedimensionaal (2D) probleem door Lagrange-vermenigvuldiger. Primaire variabelen op de subdomeinhoeken worden geëxtraheerd om de convergentiesnelheid van de iteratieve oplossing van het duale probleem te versnellen [23-26]. Een lage-rang sparsificatie benadering wordt ontwikkeld om de prestaties van de FETI-DP te verbeteren. Het gebruikt data-sparse algoritmen om de efficiëntie te verbeteren voor het oplossen van de subdomeinproblemen en het duale probleem [27, 28]. De voorgestelde methode biedt volledig gedecoupleerde subdomeinen, wat de parallelle simulatie van optische kracht voor het vangen van nanodeeltjes mogelijk maakt. De Maxwell-spanningstensor (MST) die de relatie tussen het elektromagnetische veld en de mechanische impuls onthult, wordt geadopteerd om de optische kracht te evalueren [29]. Op basis van de verkregen optische kracht kan het optische potentieel verder worden berekend voor de analyse van een stabiele vanging. Vergeleken met IEMs is de voorgestelde methode krachtiger in het omgaan met samengestelde materialen en het oplossen van het nabijveld voor SP-gebaseerde optische vanging. Vergeleken met FDTD kan de voorgestelde methode dispersief metaalmateriaal in SP-gebaseerde optische vangsystemen nauwkeurig hanteren en de trap-fout elimineren voor objecten met gekromde grens. Vergeleken met FEM is de voorgestelde methode geschikt voor grootschalige berekening van optische vanging. Verschillende voorbeelden worden geanalyseerd en numerieke resultaten demonstreren de nauwkeurigheid en efficiëntie van de voorgestelde methode voor de voorspelling en analyse van optische vanging op nanoschaal.

Methoden

FETI-DP Formuleringen

Voor de FETI-DP implementatie wordt het oorspronkelijke computationele domein Ω eerst gescheurd in een reeks non-overlappende subdomeinen Ωⁱ (i = 1, 2, 3⋯, N_s), zoals getoond in Fig. 1. In elk subdomein Ωⁱ kan een subdomein finite element systeem worden afgeleid uit de vector-golfvergelijking als

$$ \nabla \times {\mu}_r^{-1}\nabla \times {\mathbf{E}}^i-{k}_0^2{\varepsilon}_r{\mathbf{E}}^i={jk}_0{\eta}_0{\mathbf{J}}_{\mathrm{imp}}^i\kern1.08em \mathrm{in}\kern0.24em {\Omega}^i $$ (1) $$ \hat{n}\times \nabla \times {\mathbf{E}}^i+{jk}_0\hat{n}\times \left(\hat{n}\times {\mathbf{E}}^i\right)=0\kern0.96em \mathrm{on}\kern0.24em {\Gamma}_{\mathrm{ABC}}^i $$ (2)

Een domeindeelingsschema met non-overlappende subdomeinen in de FETI-DP methode. a Oorspronkelijk domein. b Verdeelde subdomeinen en gediscretiseerde meshes

waarbij Eⁱ de onbekende elektrische veld aanduidt die opgelost moet worden in \( {\Omega}^i \), k₀ en η₀ zijn de vrije-ruimte golffgetal en intrinsieke impedantie, respectievelijk, en \( {\mathbf{J}}_i^{\mathrm{imp}} \) is de geïmpresseerde stroom. \( {\Gamma}_{\mathrm{ABC}}^i \) betekent de absorberende randvoorwaarde (ABC) om de oneindige open regio af te knopen. Het moet worden opgemerkt dat k₀ vervangen moet worden door de golffimpedantie in medium als de omringende medium geen vrije-ruimte is, wat veelvoorkomend is voor optische vanging. Op de subdomeininterface Γⁱ is een aangenomen randvoorwaarde vereist om een compleet randwaardeprobleem in Ωⁱ te genereren. Hier wordt een Robin-type transmissievoorwaarde met een onbekende hulpvariabele Λ_i opgelegd als

$$ {\hat{n}}^i\times \left({\mu}_r^{-1}\nabla \times {\mathbf{E}}^i\right)+{\alpha}^i{\hat{n}}^i\times \left({\hat{n}}^i\times {\mathbf{E}}^i\right)={\boldsymbol{\Lambda}}^i\kern1.2em \mathrm{on}\kern0.36em {\Gamma}^i $$ (3)

waarbij \( {\hat{n}}^i \) de eenheidsnormale uitwaartse vector aanduidt op de subdomeininterface Γⁱ, en αⁱ een complexe parameter is die vaak kan worden gekozen als jk₀. Alle subdomeinen worden vervolgens gediscretiseerd door tetrahedrale elementen. In elk element breiden we E uit met vector basisfuncties N en onbekende elektrisch veldcoëfficiënt E als

$$ \mathbf{E}=\sum \limits_{p=1}^s{E}_p{\mathbf{N}}_p $$ (4)

waarbij s het aantal vector basisfuncties in elk tetrahedraal element aanduidt. s wordt gekozen als 6 voor traditionele lage-orde basisfunctie gebaseerd op de rand, terwijl het groter is voor hoge-orde vector basisfunctie, aangezien additionele basisfuncties gebaseerd op vlak of volume worden geïntroduceerd.

Toepassen van Galerkin's methode, kan de FEM matrixvergelijking in Ωⁱ over de onbekende elektrisch veldcoëfficiënt Eⁱ worden verkregen als

$$ \left(\begin{array}{cc}{\mathbf{K}}_{rr}^i& {\mathbf{K}}_{rc}^i\\ {}{\mathbf{K}}_{cr}^i& {\mathbf{K}}_{cc}^i\end{array}\right)\left(\begin{array}{l}{E}_r^i\\ {}{E}_c^i\end{array}\right)=\left(\begin{array}{l}{f}_r^i-{\mathbf{B}}_r^{i^T}{\lambda}^i\\ {}{f}_c^i\end{array}\right) $$ (5)

waarbij de subscript notaties c en r de hoek graad van vrijheid (DOFs) en de overige DOFs onderscheiden, wat de hoek DOFs extraheert als een primair variabele om het dual-primal (DP) schema te construeren. Hier is K de FEM systeemmatrix en f de excitatievector. B is een Booleaanse matrix die de interface DOFs van een subdomein extraheert. λ is een duale variabele gegenereerd uit het uitbreiden van Λ_i, wat ook de Lagrange-vermenigvuldiger wordt genoemd.

Dan kunnen de aangrenzende subdomeinen worden interconnected door het afdwingen van tangentiële elektrische veld en magnetische veld continuïteit op de interfaces. We assembleren alle subdomeininterfaces en elimineren alle subdomein interne onbekenden Eⁱ. Een gereduceerde globale interfacevergelijking over de duale variabele λ kan worden verkregen als

$$ \left[{\tilde{\mathbf{K}}}_{rr}+{\tilde{\mathbf{K}}}_{rc}{\tilde{\mathbf{K}}}_{cc}^{-1}{\tilde{\mathbf{K}}}_{cr}\right]\lambda ={\tilde{f}}_r-{\tilde{\mathbf{K}}}_{rc}{\tilde{\mathbf{K}}}_{cc}^{-1}{\tilde{f}}_c $$ (6)

Vergelijking (6) kan worden opgelost door iteratieve methoden, zoals de gegeneraliseerde minimale residu (GMRES) methode. \( {\tilde{\mathbf{K}}}_{cc} \) is het globale hoeksysteem, wat de iteratieve convergentie in primale ruimte kan versnellen. Geschikte preconditioner kan worden gebruikt om de iteratieve convergentiesnelheid te verbeteren, zoals benaderde inverse en onvolledige LU-decompositie. Eenmaal de duale variabele λ opgelost is, kan het elektrische veld binnen elk subdomein onafhankelijk worden geëvalueerd door (5). Voor de constructie van de globale matrix \( {\tilde{\mathbf{K}}}_{rr} \), moet men de subdomein numerieke Green's functie \( {\mathbf{Z}}_{rr}^i \) construeren met een vorm van

$$ {\mathbf{Z}}_{rr}^i={\mathbf{B}}_r^i{\left({\mathbf{K}}_{rr}^i\right)}^{-1}{{\mathbf{B}}_r^i}^T $$ (7)

waarbij de inverse van de subdomein FEM matrix \( {\left({\mathbf{K}}_{rr}^i\right)}^{-1} \) is inbegrepen. Bovendien, voor matrices \( {\tilde{\mathbf{K}}}_{rc} \), \( {\tilde{\mathbf{K}}}_{cr} \), en \( {\tilde{\mathbf{K}}}_{cc} \) en vectoren \( {\tilde{f}}_r \) en \( {\tilde{f}}_c \), is de \( {\left({\mathbf{K}}_{rr}^i\right)}^{-1} \) vereist om berekend te worden. De constructies van \( {\left({\mathbf{K}}_{rr}^i\right)}^{-1} \) in de pre-processing fase alsook hun matrix-vector producten (MVPs) in de iteratieve oplossingsfase zijn computationeel duur. Hoewel \( {\mathbf{K}}_{rr}^i \) sparse is, zijn \( {\left({\mathbf{K}}_{rr}^i\right)}^{-1} \) dicht, wat hoge computationele kosten betekent. Vervolgens wordt een lage-rang sparsificatie methode geïntroduceerd om de berekening van \( {\left({\mathbf{K}}_{rr}^i\right)}^{-1} \) te versnellen. Aangezien sommige sub-matrices in het globale interfacesysteem in lage-rang matrixvorm kunnen worden voorgesteld, kan hun berekening worden uitgevoerd met lage-rang algoritme, wat de prestaties van de FETI-DP verbetert. Het kan worden gezien dat de FETI-DP onafhankelijke subdomeinoperaties biedt, zodat het parallelle berekening kan exploiteren om de efficiëntie te verbeteren. Voor een efficiënt parallel schema is een principe van domeindeeling om het aantal DOFs in elk subdomein zo gebalanceerd mogelijk te maken. Daarom moet de grootte van subdomeinen verband houden met de mesh-discretisatiedichtheid. Meestal worden kleine subdomeinen geadopteerd in fijn gegedeelde gebieden, terwijl grote subdomeinen worden geadopteerd in grof gegedeelde gebieden.

Lage-rang Sparsificatie

Een lage-rang sparsificatie benadering wordt voorgesteld om een data-sparse manier te bieden om de FETI-DP efficiëntie te verbeteren. Hier betekent data-sparse dat deze matrices eigenlijk niet sparse zijn maar ze zijn sparse in de zin dat bepaalde sub-blokken van hen kunnen worden voorgesteld door lage-rang decompositie matrixvormen als

$$ \mathbf{M}={\mathbf{XY}}^{\mathrm{T}}\kern0.72em \left(\mathbf{M}\in {\mathrm{\mathbb{C}}}^{m\times n},\mathbf{X}\in {\mathrm{\mathbb{C}}}^{m\times k},\mathbf{Y}\in {\mathrm{\mathbb{C}}}^{n\times k}\right) $$ (8)

waarbij X en Y in volle matrixvormen zijn, en de rang k veel kleiner is dan m en n. De matrix \( {\left({\mathbf{K}}_{rr}^i\right)}^{-1} \) kan worden voorgesteld door data-sparse matrixvormen aangezien het de matrixeigenschap van een integraaloperator bezit. Dus, mits \( {\left({\mathbf{K}}_{rr}^i\right)}^{-1} \) lage-rang eigenschap bezit in gegeven subdomein, kan het efficiënt worden berekend en opgeslagen in data-sparse vormen met de lage-rang sparsificatie benadering, wat de MVPs in de iteratieve oplossing versnelt.

De processen van de lage-rang sparsificatie benadering kunnen worden verdeeld in de volgende stappen: (1) construeren een clusterboom door onderverdeling van de basisfunctieset in elk subdomein, (2) construeren een blok-clusterboom door interactie van twee clusterbomen, (3) genereren een data-sparse vorm van \( {\mathbf{K}}_{rr}^i \) door een toelaatbaarheidsvoorwaarde, (4) uitvoeren lage-rang geformatteerde algoritmen om de data-sparse representatie van \( {\left({\mathbf{K}}_{rr}^i\right)}_{\mathrm{DS}}^{-1} \) te krijgen, en (5) binnenkomen de oplossing van FETI-DP systemen door data-sparse algoritme. Geschikte preconditioner kan worden gebruikt om de oplossing te versnellen. Het moet worden opgemerkt dat de data-sparse LU-factorisatie \( {\left({\mathbf{K}}_{rr}^i\right)}_{\mathrm{DS}}={\left({\mathbf{L}}_{rr}^i\right)}_{\mathrm{DS}}{\left({\mathbf{U}}_{rr}^i\right)}_{\mathrm{DS}} \) wordt geadopteerd om de matrixinversie \( {\left({\mathbf{K}}_{rr}^i\right)}_{\mathrm{DS}}^{-1} \) te vervangen. Een geneste dissectie techniek wordt gebruikt om de efficiëntie van de lage-rang sparsificatie verder te verbeteren. De geneste dissectie gebruikt scheiders om grote off-diagonale nul sub-blokken te produceren, wat nul zal blijven tijdens de LU-factorisatie zodat het de fill-ins significant kan verminderen.

Om \( {\left({\mathbf{K}}_{rr}^i\right)}_{\mathrm{DS}} \) te genereren, construeren we eerst een clusterboom T_I door recursieve onderverdeling van de subdomein rand-gebaseerde basisfunctieset I = {1,2,……N} met behulp van bounding box. Met de geneste dissectie wordt een cluster t binnen de corresponderende bounding box verdeeld in drie opvolgers {s₁, s_sep, s₂}, waar s₁ en s₂ de indexsets zijn van de twee ontkoppelde bounding boxes en s_sep de indexset is van de scheider. Figuur 2 a toont een eenvoudig voorbeeld van dit proces. Dan kan een blok-clusterboom T_{I × I} worden geconstrueerd door interactie van twee clusterboom T_I, zoals getoond in Fig. 2 b, wat kan worden gekozen als de clusterboom van de oorspronkelijke rand-gebaseerde basisfunctieset en die van de test basisfunctieset in Galerkin's methode. Vervolgens moeten we een toelaatbaarheidsvoorwaarde introduceren gebaseerd op de geneste dissectie om volle blokken, lage-rang decompositie blokken en off-diagonale nulblokken in T_{I × I} te onderscheiden [23]. Dus kan \( {\left({\mathbf{K}}_{rr}^i\right)}_{\mathrm{DS}} \) worden geproduceerd door de corresponderende blokken te vullen met de niet-nul entries van \( {\mathbf{K}}_{rr}^i \). Tenslotte kan de data-sparse LU-factorisatie van \( {\left({\mathbf{K}}_{rr}^i\right)}_{\mathrm{DS}}={\left({\mathbf{L}}_{rr}^i\right)}_{\mathrm{DS}}{\left({\mathbf{U}}_{rr}^i\right)}_{\mathrm{DS}} \) recursief worden berekend vanaf

$$ {\mathbf{K}}_{rr}^i=\left[\begin{array}{ccc}{\mathbf{K}}_{11}& & {\mathbf{K}}_{13}\\ {}& {\mathbf{K}}_{22}& {\mathbf{K}}_{23}\\ {}{\mathbf{K}}_{31}& {\mathbf{K}}_{32}& {\mathbf{K}}_{33}\end{array}\right]=\left[\begin{array}{ccc}{\mathbf{L}}_{11}& & \\ {}& {\mathbf{L}}_{22}& \\ {}{\mathbf{L}}_{31}& {\mathbf{L}}_{32}& {\mathbf{L}}_{33}\end{array}\right]\left[\begin{array}{ccc}{\mathbf{U}}_{11}& & {\mathbf{U}}_{13}\\ {}& {\mathbf{U}}_{22}& {\mathbf{U}}_{23}\\ {}& & {\mathbf{U}}_{33}\end{array}\right] $$ (9)

Constructies van een clusterboom en een blok-clusterboom van 4 niveaus gebaseerd op geneste dissectie. a Constructie van een clusterboom door recursieve onderverdeling van rand-gebaseerde basisfunctieset I = {1,2,…18}. b Constructie van een blok-clusterboom waar witte blokken nulmatrices zijn en groene blokken volle matrices of lage-rang decompositie matrices kunnen zijn

waarbij conventionele volle matrixrekenkunde wordt vervangen door hun data-sparse tegenhangers [28]. Een adaptieve afkapingsfout ε_t wordt gebruikt om de nauwkeurigheid van lage-rang benaderingen te controleren. De verkregen LU-factoren \( {\left({\mathbf{L}}_{rr}^i\right)}_{\mathrm{DS}} \) en \( {\left({\mathbf{U}}_{rr}^i\right)}_{\mathrm{DS}} \) worden opgeslagen en gebruikt om \( {\mathbf{Z}}_{rr}^i \) te construeren door

$$ {\mathbf{Z}}_{rr}^i={\mathbf{B}}_r^i{\left({\mathbf{U}}_{rr}^i\right)}_{\mathrm{DS}}^{-1}{\left({\mathbf{L}}_{rr}^i\right)}_{\mathrm{DS}}^{-1}{\mathbf{B}}_r^i $$ (10)

waarbij \( {\mathbf{B}}_r^i{\left({\mathbf{U}}_{rr}^i\right)}_{\mathrm{DS}}^{-1} \) en \( {\left({\mathbf{L}}_{rr}^i\right)}_{\mathrm{DS}}^{-1}{\mathbf{B}}_r^i \) kunnen worden berekend door data-sparse boven- en onderdriehoekige solver. De \( {\left({\mathbf{L}}_{rr}^i\right)}_{\mathrm{DS}} \), \( {\left({\mathbf{U}}_{rr}^i\right)}_{\mathrm{DS}} \), en \( {\mathbf{Z}}_{rr}^i \) gaan in de FETI-DP berekening met data-sparse voorwaartse en achterwaartse substituties (FBSs) en data-sparse MVPs.

Optische Kracht en Potentieel

Volgens elektrodynamicatheorie kan de optische kracht worden geëvalueerd door de Maxwell-spanningstensor (MST) die de relatie tussen elektromagnetisch veld en mechanische impuls onthult [29]. Eenmaal de elektromagnetische veldverdeling rond het object verkregen is, kan de optische kracht worden berekend door MST te integreren over een gesloten oppervlak dat het object omringt. Op basis van de verkregen elektrisch veldverdeling kan de MST op elke coördinaat worden geconstrueerd door

$$ \overleftrightarrow{\mathbf{T}}=\frac{1}{2}\operatorname{Re}\left[\varepsilon {\mathbf{EE}}^{\ast }+\mu {\mathbf{HH}}^{\ast }-\frac{1}{2}\left(\varepsilon {\left|\mathbf{E}\right|}^2+\mu {\left|\mathbf{H}\right|}^2\right)\overleftrightarrow{\mathbf{I}}\right] $$ (11)

waarbij de superscript asterisk de geconjugeerde van elektrisch veld of magnetisch veld aanduidt, ε zijn μ de permittiviteit en permeabiliteit, en \( \overleftrightarrow{\mathbf{I}} \) is een 3 × 3 identiteitsmatrix. Door het outer product van vectoren kan de tensorvorm van \( \overleftrightarrow{\mathbf{T}} \) worden geschreven als

$$ \overleftrightarrow{\mathrm{T}}=\left[\begin{array}{lll}{T}_{xx}& {T}_{xy}& {T}_{xz}\\ {}{T}_{yx}& {T}_{yy}& {T}_{yz}\\ {}{T}_{zx}& {T}_{zy}& {T}_{zz}\end{array}\right]=\left[\begin{array}{ccc}\varepsilon {E}_x{E}_x^{\ast }+\mu {H}_x{H}_x^{\ast }-\frac{\varepsilon {\left|\mathbf{E}\right|}^2+\mu {\left|\mathbf{H}\right|}^2}{2}& \varepsilon {E}_x{E}_y^{\ast }+\mu {H}_x{H}_y^{\ast }& \varepsilon {E}_x{E}_z^{\ast }+\mu {H}_x{H}_z^{\ast}\\ {}\varepsilon {E}_y{E}_x^{\ast }+\mu {H}_y{H}_x^{\ast }& \varepsilon {E}_y{E}_y^{\ast }+\mu {H}_y{H}_y^{\ast }-\frac{\varepsilon {\left|\mathbf{E}\right|}^2+\mu {\left|\mathbf{H}\right|}^2}{2}& \varepsilon {E}_y{E}_z^{\ast }+\mu {H}_y{H}_z^{\ast}\\ {}\varepsilon {E}_z{E}_x^{\ast }+\mu {H}_z{H}_x^{\ast }& \varepsilon {E}_z{E}_y^{\ast }+\mu {H}_z{H}_y^{\ast }& \varepsilon {E}_z{E}_z^{\ast }+\mu {H}_z{H}_z^{\ast }-\frac{\varepsilon {\left|\mathbf{E}\right|}^2+\mu {\left|\mathbf{H}\right|}^2}{2}\end{array}\right] $$ (12)

waarbij de subscript x, y, z de componenten in drie richtingen aanduidt. Volgens de uitbreiding van E beschreven in (4), kunnen de entries van MST T_mn (m, n = x, y, z) worden geconverteerd naar uitgebreide vormen in de FETI-DP berekening als

$$ {\displaystyle \begin{array}{l}{T}_{mn}=\sum \limits_{p,q=1}^s{E}_p{E}_q\left\{\varepsilon {\left({\mathbf{N}}_p\right)}_m{\left({\mathbf{N}}_q^{\ast}\right)}_n-\frac{1}{\omega^2\mu }{\left(\nabla \times {\mathbf{N}}_p\right)}_m{\left(\nabla \times {\mathbf{N}}_q^{\ast}\right)}_n\right.\\ {}\kern1.75em \left.-\frac{1}{2}\left[\varepsilon \left({\mathbf{N}}_p\right)\left({\mathbf{N}}_q^{\ast}\right)-\frac{1}{\omega^2\mu}\left(\nabla \times {\mathbf{N}}_p\right)\left(\nabla \times {\mathbf{N}}_q^{\ast}\right)\right]\right\}\kern1.75em \mathrm{if}\ m=n.\end{array}} $$ (13) $$ {T}_{mn}=\sum \limits_{p,q=1}^s{E}_p{E}_q\left\{\varepsilon {\left({\mathbf{N}}_p\right)}_m{\left({\mathbf{N}}_q^{\ast}\right)}_n-\frac{1}{\omega^2\mu }{\left(\nabla \times {\mathbf{N}}_p\right)}_m{\left(\nabla \times {\mathbf{N}}_q^{\ast}\right)}_n\right\}\kern1.25em \mathrm{if}\ m\ne n. $$ (14)

waarbij ω de hoekfrequentie is; N en s zijn beschreven in Verg. (4).

Tenslotte kan de optische kracht uitgeoefend op het object worden berekend door de MST te integreren over elk gesloten oppervlak dat het omringt door

$$ \mathbf{F}={\oint}_S\left(\overleftrightarrow{\mathbf{T}}\cdot \hat{n}\right)\ dS. $$ (15)

Let op dat de berekening van optische kracht ook in parallel geïmplementeerd kan worden, aangezien de integraal van de MST wordt toegewezen aan corresponderende subdomeinen. Voor een stabiele optische vanging is een van de hoofdvoorwaarden dat de gradientkracht groter moet zijn dan de verstrooiingskracht. Met andere woorden, de richting van de totale kracht moet identiek zijn aan die van de gradientkracht, wat altijd wijst naar de positie waar de elektrische veldintensiteit het sterkst is.

Het optische potentieel is een andere aantrekkelijke parameter die de stabiliteit van de optische vanging onthult. Op basis van de verkregen optische kracht kan de optische potentiealdiepte U op positie r₀ worden berekend door

$$ \mathbf{U}\left({r}_0\right)=-{\int}_{\infty}^{r_0}\mathbf{F}\left(\mathbf{r}\right)\cdot \mathbf{r}, $$ (16)

waarbij de subscript ∞ oneindigheid aanduidt gedefinieerd als het referentiepunt met nul potentieel. De waarde van U kan worden voorgesteld door k_BT, waar k_B de Boltzmannconstante aanduidt van 1,3806488 × 10⁻²³ J/K en T de omgevingstemperatuur is. In het algemeen kan het deeltje de Brownse beweging in oplossing overwinnen en stabiel gevangen worden wanneer U > 1 k_BT voldaan is. Anders kan het deeltje niet stabiel gevangen worden. Aangezien de totale optische kracht de conservatieve gradientkrachtcomponent en de niet-conservatieve verstrooiingskrachtcomponent omvat, is de totale optische kracht F uit (15) niet-conservatief [30, 31]. Echter, mits de beweging van het nanodeeltje beperkt is tot één dimensie, levert dit een eenduidige definitie van een optisch potentieel op uit (16), zelfs al is de totale optische kracht niet-conservatief.

Resultaten en Discussie

Drie voorbeelden worden gepresenteerd om de effectiviteit van de voorgestelde methode te demonstreren. Aangezien edele metalen vaak worden gebruikt om de oppervlakteplasmon op te wekken, selecteren we representatieve goud en zilver materialen voor de analyses. Het eerste voorbeeld berekent de optische kracht van zilver nanodeeltje om de nauwkeurigheid van de voorgestelde FETI-DP methode te verifiëren. Het tweede en derde voorbeeld simuleren en bespreken de optische vanging van goud nanodeeltjes. Voor alle voorbeelden wordt het oneindige domein afgekapt met ABC, en de afstanden tussen de ABC en de objecten worden ingesteld op één golfengte, wat voldoende is om een acceptabele nauwkeurigheid te bereiken. Alle berekeningen worden uitgevoerd op een Dell workstation uitgerust met 3,6 GHz Intel Xeon processoren.

Zilver Nanolaksel

Een zilver nanolaksel object wordt eerst overwogen om de nauwkeurigheid en efficiëntie van de voorgestelde FETI-DP methode in het voorspellen van optische kracht te testen. Figuur 3 a en b presenteert de configuratie en afmetingen. De constitutieve parameters van zilver zijn alle gemeten waarden afkomstig uit [32]. Om het FETI-DP schema te implementeren, wordt het hele analysedomein eerst verdeeld in 24 subdomeinen. Dichtere meshes zijn vereist nabij het metaaloppervlak om het plasmonische lokale veldversterkingseffect te modelleren. Tetrahedrale elementen worden geadopteerd voor de discretisatie, wat leidt tot in totaal 6,9 × 10⁵ onbekenden, inclusief 4,1 × 10⁴ duale onbekenden en 313 hoekonbekenden. Het invallicht belicht langs de richting van +z, terwijl de polarisatierichting van elektrisch veld −x is.

Configuratie van een zilver nanolaksel structuur. a 3D-weergave. b Frontweergave en afmetingen, waar R = 30 nm en h = 60 nm

Eerst veranderen we de golfengte van invallicht λ van 200 nm tot 400 nm om de optische krachten uitgeoefend op de nanolaksel te simuleren. Aangezien de FETI-DP werkt in frequentiedomein, worden de optische krachten berekend bij 15 sampling frequentiepunt. Figuur 4 toont de berekende curve van optische krachten uitgeoefend op de zilver nanolaksel. Om de nauwkeurigheid van de FETI-DP aan te duiden, worden de optische krachtresultaten van de FETI-DP vergeleken met die van de commerciële software Lumerical FDTD Solutions [33], en goede overeenstemming kan worden waargenomen.

Resultaten van de optische krachten uitgeoefend op de zilver nanolaksel, varierend met de golfengte λ van invallicht, inclusief de resultaten van de FETI-DP en de commerciële software FDTD oplossingen

Dan wordt de prestatie van FETI-DP getest voor verschillende aantallen subdomeinen. We verhogen het aantal subdomeinen van 4 tot 24 door de discretisatiedichtheid te behouden. We wijzen elke processor toe om met één subdomein te dealen. Tabel 1 rapporteert de tijd gebruikt voor de constructie van globale interface Verg. (6) en de totale oplossingstijd. Het kan worden gezien dat de FETI-DP parallelle computerbronnen volledig kan benutten en de oplossingsefficiëntie significant kan verbeteren. Bovendien wordt de nauwkeurigheid van de FETI-DP met het toenemende aantal subdomeinen ook onderzocht en gerapporteerd in Tabel 1. Hier wordt de nauwkeurheid gedefinieerd door de 2-norm relatieve fout van de optische kracht als δ_OF = ‖OF_i − OF_ref‖/‖OF_ref‖, waar OF_i de optische kracht is met i subdomeinen en OF_ref de referentie optische kracht aanduidt met twee subdomeinen. Het kan worden gezien dat de nauwkeurigheid bijna constant blijft met het toenemende aantal subdomeinen.

Goud Nanospheren Dimer

Het tweede voorbeeld analyseert de optische vanging van een goud nanospher door gebruik te maken van een goud nanospheren dimer. De plasmonische effecten bij de dimerkloof kunnen de optische kracht effectief versterken voor het vangen van nanodeeltjes. Figuur 5 a en b geeft de configuratie en afmetingen van dit systeem. De constitutieve parameters van goud zijn alle gemeten waarden afkomstig uit [32]. Het omringende medium is water met een relatieve brekingsindex van n = 1,33. Het invallicht is een vlakke golf met een vermogen van 10 mW/μm², de elektrisch veld polarisatierichting is +x, en de invalsrichting is −z. De optische kracht uitgeoefend op het object nanospher wordt berekend door de FETI-DP methode. Voor de FETI-DP implementatie wordt het hele computationele domein verdeeld in 32 subdomeinen en gediscretiseerd door tetrahedrale meshes, wat resulteert in 3,5 × 10⁶ onbekenden, inclusief 1,6 × 10⁵ duale onbekenden en 1738 hoekonbekenden.

Configuratie van een optisch vangingssysteem van een goud nanospheren dimer in water. a 3D-weergave. b Frontweergave en afmetingen, waar R = 25 nm, r = 5 nm, en g = 2 nm

Eerst testen we de parallelle prestaties van de voorgestelde FETI-DP door gebruik te maken van verschillende aantallen processoren. Tabel 2 rapporteert de oplossingstijd voor Verg. (6) alsook de totale oplossingstijd. Bovendien worden de versnellingen voor de parallelle berekening ook in Tabel 2 gegeven. Hier wordt de versnelling gedefinieerd door

$$ \mathrm{Speed}\ \mathrm{up}=\raisebox{1ex}{${T}_1$}\!\left/ \!\raisebox{-1ex}{${T}_{N_p}$}\right. $$ (17)

waarbij \( {T}_{N_p} \) de totale kloktijd aanduidt met N_p processoren. Het kan worden gezien dat de FETI-DP de oplossingsefficiëntie significant verbetert en goede parallelle versnelling vertoont. Voor dit grote aantal onbekenden is het totale geheugengebruik van alle processoren slechts 57,2 GB.

Dan wordt de effectiviteit van de lage-rang sparsificatie benadering onderzocht. Met de lage-rang sparsificatie kan de subdomeinmatrix worden gefactoriseerd door data-sparse algoritme en opgeslagen als data-sparse matrices. De constructietijd en geheugengebruik zijn slechts 18 s en 0,5 GB, terwijl ze 67 s en 1,7 GB zijn door conventioneel matrixalgoritme. Het kan worden gezien dat we 72% tijdsbesparing en 70% geheugencompressie krijgen. Gerelateerd aan het geheugengebruik kunnen de daaropvolgende MVPs ook 70% tijdsbesparing krijgen.

Vervolgens wordt de FETI-DP getest voor de optische krachtberekening met de golfengte λ varierend van 277 nm tot 818 nm. In de praktijk zijn analyses van optische kracht onder invallicht van verschillende golfengtes vaak noodzakelijk voor het zoeken naar de plasmonische resonantiegolfengte, waar drastische veldversterking optreedt en de sterkste optische kracht kan worden verkregen. Twee gevallen worden overwogen met de nanospher gelegen op (0, 0, 20 nm) en (0, 0, − 20 nm). Figuur 6 a en b plot de berekende optische krachten uitgeoefend op de nanospher voor verschillende λ. Het kan worden gezien dat de maximum optische kracht optreedt bij λ = 472 nm, wat de plasmonische resonantiegolfengte is. De optische kracht bij deze resonantiegolfengte versterkt met bijna 40 keer ten opzichte van die bij niet-resonantiegolfengte. Bovendien wijst de optische kracht altijd naar de dimerkloof, zoals getoond in Fig. 6, waar de elektrische veldintensiteit het sterkst is. Het is ook de richting van gradientkracht om het object te vangen. Figuur 7 a en b toont de berekende elektrische veldversterkingsverdelingen bij de niet-resonantiegolfengte van λ = 300 nm en de resonantiegolfengte van λ = 472 nm, respectievelijk. Het kan worden gezien dat de elektrische veldintensiteit is verhoogd met bijna 250 keer door het plasmonische resonantieeffect.

Berekende resultaten van optische krachten uitgeoefend op de nanospher in het systeem van goud nanospheren dimer, varierend met de golfengte λ van invallicht. a Het object nanospher is gelegen op (0, 0, 20 nm). b Het object nanospher is gelegen op (0, 0, − 20 nm)

De elektrische veldversterkingsverdelingen op het xoz vlak voor het systeem van goud nanospheren dimer. a λ = 300 nm (niet-resonantiegolfengte). b λ = 472 nm (resonantiegolfengte)

Bovendien worden de optische kracht en optisch potentieel berekend met de nanospher bewegend van (0, 0, − 30 nm) tot (0, 0, − 17 nm) langs de z-as. Aangezien de meest typische en interessante gedragingen van vangkrachten en potentiellen die zijn die werkend langs z-richting, overwegen we hier de axiële vangpotentieel door integratie langs de z-as. Omdat de beweging van het nanodeeltje beperkt is tot één dimensie, is de definitie van een optisch potentieel eenduidig uit (16), zelfs al is de totale optische kracht uit (15) niet-conservatief. Zoals getoond in Fig. 8 a, b, met de nanospher bewegend naar de dimerkloof, nemen de optische kracht en optische potentiealdiepte duidelijk toe. Op de positie van (0, 0, − 17 nm) wordt een optische potentiealdiepte van 4,6 k_BT geproduceerd, wat voldoende is om de Brownse beweging in water te overwinnen om stabiele optische vanging te bereiken.

De optische krachten en optische potentiëlen uitgeoefend op de nanospher in het systeem van goud nanospheren dimer, wanneer de nanospher beweegt van (0, 0, − 30 nm) tot (0, 0, − 17 nm). a De optische krachten. b De optische potentiëlen

Tenslotte testen we de effecten van het dielektrische substraat voor dit voorbeeld. De optische krachten worden berekend met en zonder substraat, respectievelijk. Voor beide twee gevallen is de nanospher gelegen op (0, 0, − 20 nm) en de invalgolfengte wordt gekozen als de resonantiegolfengte. Voor het geval zonder substraat is het berekende resultaat van de optische kracht |F₀| = 0,769 pN. Voor het geval met een substraat wordt de goud nanospheren dimer op een dielektrisch substraat geplaatst met een dikte van 60 nm en een relatieve permittiviteit van ε_r = 2,25. Het berekende resultaat van de optische kracht is |F₁| = 0,761 pN. De relatieve fout tussen deze twee resultaten van optische krachten is ongeveer 1,0 × 10⁻², wat gedefinieerd is als |F₁ − F₀|/|F₀|.

Goud Afgeknotte Kegel Dimer

Het derde voorbeeld gaat over de optische vanging van een goud nanospher door gebruik te maken van een goud afgeknotte kegel dimer. Figuur 9 geeft de configuratie en afmetingen van dit systeem. De constitutieve parameters van goud zijn afkomstig uit [32]. Het dielektrische substraat heeft een relatieve permittiviteit van ε_r = 2,25. Het omringende medium is water met een relatieve brekingsindex van n = 1,33. Het invallicht is vlakke golf met een vermogen van 10 mW/μm², de elektrisch veld polarisatierichting is +x, en de invalsrichting is −z. Het hele computationele domein wordt verdeeld in 32 subdomeinen en gediscretiseerd door tetrahedrale meshes, wat leidt tot 3,1 × 10⁶ onbekenden, inclusief 1,3 × 10⁵ duale onbekenden en 1227 hoekonbekenden.

Configuratie van een optisch vangingssysteem van een goud afgeknotte kegel dimer gebaseerd op een dielektrisch substraat in water. a 3D-weergave. b Frontweergave en afmetingen, waar UR = 20 nm, LR = 30 nm, h = 35 nm, en g = 2 nm

Eerst analyseren we de optische krachten door λ te veranderen van 277 nm tot 818 nm. Figuur 10 plot de berekende optische krachten uitgeoefend op de nanospher voor verschillende λ. De nanospher is gelegen op (0, 0, 35 nm). Het kan worden gezien dat de maximum optische kracht optreedt bij λ = 464 nm, wat de plasmonische resonantiegolfengte is, en de optische kracht hier wordt versterkt met bijna 30 keer bij niet-resonantiegolfengte. Bovendien wijst de totale optische kracht altijd naar −z, zoals getoond in Fig. 10, wat de richting is van de gradientkracht. Dit bevestigt dat de gradientkracht groter is dan de verstrooiingskracht, wat een van de voorwaarden is dat de nanospher stabiel gevangen kan worden. Figuur 11 a en b presenteert de berekende elektrische veldverdelingen bij de niet-resonantiegolfengte van λ=300 nm en de resonantiegolfengte van λ = 464 nm, respectievelijk. Het kan worden gezien dat elektrische veldintensiteit is verhoogd met bijna 500 keer door de gelokaliseerde oppervlakteplasmonresonantie.

Berekende resultaten van optische krachten uitgeoefend op de nanospher in het systeem van goud afgeknotte kegel dimer, varierend met λ. De nanospher is gelegen op (0, 0, 35 nm)

De elektrische veldversterkingsverdelingen op het xoz vlak voor het systeem van goud afgeknotte kegel dimer. a λ= 300 nm (niet-resonantiegolfengte). b λ= 464 nm (resonantiegolfengte)

Dan wordt de locatie van de nanospher gewijzigd naar 0, 5 en 35 nm om de optische kracht te observeren. Figuur 12 geeft de berekende optische krachten uitgeoefend op de nanospher, waar duidelijke y-component van optische kracht kan worden waargenomen, terwijl grotere z-component van optische kracht bestaat. De totale optische kracht wijst nog steeds naar de positie met de sterkste elektrische veld om de nanospher te vangen.

Berekende resultaten van optische krachten uitgeoefend op de nanospher in het systeem van goud afgeknotte kegel dimer varierend λ. De nanospher is gelegen op (0, 5 nm, 35 nm)

Verder analyseren we het optische potentieel met de nanospher bewegend van (0, 0, 50 nm) tot (0, 0, 20 nm) langs de z-as. Hier overwegen we de axiële vangpotentieel langs z-richting, wat de beweging van het nanodeeltje beperkt tot één dimensie en leidt tot een eenduidige definitie van optisch potentieel. Zowel de optische kracht als potentieel worden berekend. Zoals kan worden waargenomen uit Fig. 13 a, b, met de nanospher bewegend naar de dimerkloof, nemen de optische kracht en de optische potentiealdiepte duidelijk toe. Bij (0, 0, 20 nm) wordt een optische potentiealdiepte van 3,8 k_BT verkregen, wat voldoende is om de Brownse beweging in water te overwinnen om stabiele optische vanging te bereiken.

De optische krachten en optische potentiëlen uitgeoefend op de nanospher in het systeem van goud afgeknotte kegel dimer, wanneer de nanospher beweegt van (0, 0, 50 nm) tot (0, 0, 20 nm). a De optische krachten. b De optische potentiëlen

Tenslotte testen we de computationele kosten van de FETI-DP door het aantal onbekenden te veranderen van 1,0 miljoen tot 3,2 miljoen gebaseerd op verschillende mesh-grootte. In de praktijk zijn tests onder verschillende mesh-dichtheid meestal noodzakelijk om aan verschillende nauwkeurigheidsvereisten te voldoen. Zo'n grootschalig complex probleem brengt grote uitdagingen voor conventionele numerieke methoden. Echter, de FETI-DP kan dit probleem gemakkelijk handeren. Tweeëndertig processoren worden gebruikt voor de FETI-DP simulatie, terwijl elke processor met een subdomein omgaat. Tabel 3 rapporteert de computationele kosten van de FETI-DP. Het kan worden gezien dat de FETI-DP hoge simulatie-efficiëntie en laag geheugengebruik vertoont.

Conclusie

Een FETI-DP methode gecombineerd met lage-rang sparsificatie wordt voorgesteld voor de voorspelling en analyse van optische vanging van metaal nanodeeltjes. De voorgestelde methode biedt volledig gedecoupleerde subdomeinproblemen, wat een grootschalig complex probleem converteert in een reeks kleinschalige eenvoudige problemen. Het is goed geschikt voor parallelle berekening en kan de efficiëntie van numerieke simulatie significant verbeteren. Voorbeelden demonstreren dat de voorgestelde methode uitstekende prestaties van grootschalige berekening vertoont en goed geschikt is voor de snelle en nauwkeurige simulatie van optische vanging op nanoschaal.

Beschikbaarheid van Data en Materialen

Alle data gegenereerd of geanalyseerd tijdens deze studie zijn opgenomen in dit artikel.

Afkortingen

ABC:

Absorberende randvoorwaarde

DOF:

Graden van vrijheid

FDTD:

Finite-difference time-domain

FEM:

Finite element method

FETI-DP:

Dual-primal finite element tearing and interconnecting

MST:

Maxwell-spanningstensor

MVP:

Matrix-vector product