%% График от E^2 и E^4 в зависимости от концентрации (меняем eplion) %% FDTD Lumerical ? %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% %This is the LaTeX ARTICLE template for RSC journals %Copyright The Royal Society of Chemistry 2014 %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% \documentclass[twoside,twocolumn,9pt]{article} \usepackage{extsizes} \usepackage[super,sort&compress,comma]{natbib} \usepackage[version=3]{mhchem} \usepackage[left=1.5cm, right=1.5cm, top=1.785cm, bottom=2.0cm]{geometry} \usepackage{balance} \usepackage{widetext} \usepackage{times,mathptmx} \usepackage{sectsty} \usepackage{graphicx} \usepackage{lastpage} \usepackage[format=plain,justification=raggedright,singlelinecheck=false,font={stretch=1.125,small,sf},labelfont=bf,labelsep=space]{caption} \usepackage{float} \usepackage{fancyhdr} \usepackage{fnpos} \usepackage{tabularx} \usepackage[english]{babel} \usepackage{array} \usepackage{droidsans} \usepackage{charter} \usepackage[T1]{fontenc} \usepackage[usenames,dvipsnames]{xcolor} \usepackage{setspace} \usepackage[compact]{titlesec} %%% Please don't disable any packages in the preamble, as this may %%% cause the template to display incorrectly.%%% \usepackage{amsmath} \newcommand{\red}[1]{{\color{red} #1}} \usepackage{epstopdf}%This line makes .eps figures into .pdf - please %comment out if not required. \definecolor{cream}{RGB}{222,217,201} \begin{document} \pagestyle{fancy} \thispagestyle{plain} \fancypagestyle{plain}{ %%%HEADER%%% \fancyhead[C]{\includegraphics[width=18.5cm]{head_foot/header_bar}} \fancyhead[L]{\hspace{0cm}\vspace{1.5cm}\includegraphics[height=30pt]{head_foot/journal_name}} \fancyhead[R]{\hspace{0cm}\vspace{1.7cm}\includegraphics[height=55pt]{head_foot/RSC_LOGO_CMYK}} \renewcommand{\headrulewidth}{0pt} } %%%END OF HEADER%%% %%%PAGE SETUP - Please do not change any commands within this %%% section%%% \makeFNbottom \makeatletter \renewcommand\LARGE{\@setfontsize\LARGE{15pt}{17}} \renewcommand\Large{\@setfontsize\Large{12pt}{14}} \renewcommand\large{\@setfontsize\large{10pt}{12}} \renewcommand\footnotesize{\@setfontsize\footnotesize{7pt}{10}} \makeatother \renewcommand{\thefootnote}{\fnsymbol{footnote}} \renewcommand\footnoterule{\vspace*{1pt}% \color{cream}\hrule width 3.5in height 0.4pt \color{black}\vspace*{5pt}} \setcounter{secnumdepth}{5} \makeatletter \renewcommand\@biblabel[1]{#1} \renewcommand\@makefntext[1]% {\noindent\makebox[0pt][r]{\@thefnmark\,}#1} \makeatother \renewcommand{\figurename}{\small{Fig.}~} \sectionfont{\sffamily\Large} \subsectionfont{\normalsize} \subsubsectionfont{\bf} \setstretch{1.125} %In particular, please do not alter this line. \setlength{\skip\footins}{0.8cm} \setlength{\footnotesep}{0.25cm} \setlength{\jot}{10pt} \titlespacing*{\section}{0pt}{4pt}{4pt} \titlespacing*{\subsection}{0pt}{15pt}{1pt} %%%END OF PAGE SETUP%%% %%%FOOTER%%% \fancyfoot{} \fancyfoot[LO,RE]{\vspace{-7.1pt}\includegraphics[height=9pt]{head_foot/LF}} \fancyfoot[CO]{\vspace{-7.1pt}\hspace{13.2cm}\includegraphics{head_foot/RF}} \fancyfoot[CE]{\vspace{-7.2pt}\hspace{-14.2cm}\includegraphics{head_foot/RF}} \fancyfoot[RO]{\footnotesize{\sffamily{1--\pageref{LastPage} ~\textbar \hspace{2pt}\thepage}}} \fancyfoot[LE]{\footnotesize{\sffamily{\thepage~\textbar\hspace{3.45cm} 1--\pageref{LastPage}}}} \fancyhead{} \renewcommand{\headrulewidth}{0pt} \renewcommand{\footrulewidth}{0pt} \setlength{\arrayrulewidth}{1pt} \setlength{\columnsep}{6.5mm} \setlength\bibsep{1pt} %%%END OF FOOTER%%% %%%FIGURE SETUP - please do not change any commands within this %%% section%%% \makeatletter \newlength{\figrulesep} \setlength{\figrulesep}{0.5\textfloatsep} \newcommand{\topfigrule}{\vspace*{-1pt}% \noindent{\color{cream}\rule[-\figrulesep]{\columnwidth}{1.5pt}} } \newcommand{\botfigrule}{\vspace*{-2pt}% \noindent{\color{cream}\rule[\figrulesep]{\columnwidth}{1.5pt}} } \newcommand{\dblfigrule}{\vspace*{-1pt}% \noindent{\color{cream}\rule[-\figrulesep]{\textwidth}{1.5pt}} } \makeatother %%%END OF FIGURE SETUP%%% %%%TITLE, AUTHORS AND ABSTRACT%%% \twocolumn[ \begin{@twocolumnfalse} \vspace{3cm} \sffamily \begin{tabular}{m{4.5cm} p{13.5cm} } \includegraphics{head_foot/DOI} & \noindent\LARGE{\textbf{Photogenerated Electron-Hole Plasma Induced Symmetry Breaking in Spherical Silicon Nanoparticle}} \\%Article title goes here instead of the text "This is the title" \vspace{0.3cm} & \vspace{0.3cm} \\ & \noindent\large{Anton Rudenko,$^{\ast}$\textit{$^{a}$} Konstantin Ladutenko,\textit{$^{b}$} Sergey Makarov\textit{$^{b}$} and Tatiana E. Itina\textit{$^{a}$$^{b}$} \textit{$^{a}$~Laboratoire Hubert Curien, UMR CNRS 5516, University of Lyon/UJM, 42000, Saint-Etienne, France } \textit{$^{b}$~ITMO University, Kronverksiy pr. 49, Saint-Petersburg, Russia} } \\%Author names go here instead of "Full name", etc. \includegraphics{head_foot/dates} & \noindent\normalsize {The concept of nonlinear all-dielectric nanophotonics based on high refractive index (e.g., silicon) nanoparticles supporting magnetic optical response has recently emerged as a powerful tool for ultrafast all-optical modulation at nanoscale. A strong modulation can be achieved via photo-generation of dense electron-hole plasma in the regime of simultaneous excitation of electric and magnetic Mie resonances, resulting in an effective transient reconfiguration of nanoparticle scattering properties. However, only homogeneous plasma generation was previously considered in a photo-excited nanoparticle, remaining unexplored any effects related to the plasma-induced optical inhomogeneities. Here we examine these effects by using 3D numerical modeling of coupled electrodynamic and material ionization equations. Based on the simulation results, we observed a deeply subwavelength plasma-induced nanopatterning of spherical silicon nanoparticles. In particular, we revealed strong symmetry breaking in the initially symmetrical nanoparticle, which arises during ultrafast photoexcitation near the magnetic dipole resonance. The proposed ultrafast breaking of the nanoparticle symmetry paves the way to the novel opportunities for nonlinear optical nanodevices.} \end{tabular} \end{@twocolumnfalse} \vspace{0.6cm} ]%%%END OF TITLE, AUTHORS AND ABSTRACT%%% %%%FONT SETUP - please do not change any commands within this section \renewcommand*\rmdefault{bch}\normalfont\upshape \rmfamily \section*{} \vspace{-1cm} %%%FOOTNOTES%%% \footnotetext{\textit{$^{a}$~Univ Lyon, UJM-St-Etienne, CNRS UMR 5516, F-42000, Saint-Etienne, France}} \footnotetext{\textit{$^{b}$~ITMO University, Kronverksiy pr. 49, St. Petersburg, Russia}} % Please use \dag to cite the ESI in the main text of the article. % If you article does not have ESI please remove the the \dag symbol % from the title and the footnotetext below. % \footnotetext{\dag~Electronic Supplementary Information (ESI) % available: [details of any supplementary information available % should be included here]. See DOI:10.1039/b000000x/} %additional % addresses can be cited as above using the lower-case letters, c, d, % e... If all authors are from the same address, no letter is required % \footnotetext{\ddag~Additional footnotes to the title and authors can % be included \emph{e.g.}\ `Present address:' or `These authors % contributed equally to this work' as above using the symbols: \ddag, % \textsection, and \P. Please place the appropriate symbol next to the % author's name and include a \texttt{\textbackslash footnotetext} entry % in the the correct place in the list.} %%%END OF FOOTNOTES%%% %%%MAIN TEXT%%%% \section{Introduction} All-dielectric nonlinear nanophotonics based on high refractive index dielectric materials has become prospective paradigm in modern optics, owing to recent advances in harmonics generation~\cite{shcherbakov2014enhanced, yang2015nonlinear, makarov2016self, shorokhov2016multifold, makarov2017efficient, makarov2017light} and ultrafast all-optical modulation~\cite{iyer2015reconfigurable, makarov2015tuning, shcherbakov2015ultrafast, yang2015nonlinear, baranov2016nonlinear, baranov2016tuning, shcherbakov2017ultrafast}. In fact, all-dielectric nanoantennas and metasurfaces possess much smaller parasitic Joule losses at high intensities as compared with their plasmonic counterparts, whereas their nonlinear properties are comparable. More importantly, the unique properties of nonlinear all-dielectric nanodevices are due to the existence of both electric and magnetic optical resonances in visible and near IR ranges~\cite{kuznetsov2016optically}. For instance, even slight variation of dielectric permittivity around optical resonances leads to significant changes in optical properties (transmittance or reflectance) of all-dielectric nanoantennas~\cite{makarov2015tuning, baranov2016nonlinear, baranov2016tuning} and metasurfaces~\cite{iyer2015reconfigurable, shcherbakov2015ultrafast, yang2015nonlinear, shcherbakov2017ultrafast}. \begin{figure}[t] \centering \includegraphics[width=0.75\linewidth]{Concept.pdf} \caption{Schematic illustration of electron-hole plasma 2D and 1D distributions in silicon nanoparticle around a magnetic resonance.} \label{fgr:concept} \end{figure} In previous works on all-dielectric nonlinear nanostructures, the building blocks (nanoparticles) were considered as objects with dielectric permittivity \textit{homogeneously} distributed over nanoparticle (NP). Therefore, in order to manipulate the propagation angle of the transmitted light it was proposed to use complicated nanostructures with reduced symmetry~\cite{albella2015switchable, baranov2016tuning, shibanuma2016unidirectional}. On the other hand, plasma explosion imaging technique~\cite{Hickstein2014} revealed \textit{in situ} strongly asymmetrical electron-hole plasma (EHP) distribution in various dielectric NPs during their pumping by femtosecond laser pulses. Therefore, local permittivity in the strongly photoexcited NPs can be significantly inhomogeneous, and symmetry of nanoparticles can be reduced. %The forward ejection of ions was attributed in this case to a nanolensing effect inside the NP and to intensity enhancement as low as $10\%$ on the far side of the NP. Much stronger enhancements can be achieved near electric and magnetic dipole resonances excited in single semiconductor NPs, such as silicon (Si), germanium (Ge) etc. In this Letter, we show theoretically that ultra-fast photo-excitation in a spherical silicon NP leads to a strongly inhomogeneous EHP distribution, as is shown schematically in Fig.~\ref{fgr:concept}. To reveal and analyze this effect, we perform a full-wave numerical simulation. We consider an intense femtosecond (\textit{fs}) laser pulse to interact with a silicon NP supporting Mie resonances and two-photon EHP generation. In particular, we couple finite-difference time-domain (FDTD) method used to solve three-dimensional Maxwell equations with kinetic equations describing nonlinear EHP generation. Three-dimensional transient variation of the material dielectric permittivity is calculated for NPs of several sizes. The obtained results propose a novel strategy to create complicated non-symmetrical nanostructures by using single photo-excited spherical silicon NPs. Moreover, we show that a dense EHP can be generated at deeply subwavelength scale ($< \lambda / 10$) supporting the formation of small metalized parts inside the NP. In fact, such effects transform a dielectric NP to a hybrid metall-dielectric one strongly extending functionality of the ultrafast optical nanoantennas. %Plan: %\begin{itemize} %\item Fig.3: Temporal evolution of EHP in NP with fixed diameter (at %MD) at different intensities, in order to show possible regimes of %plasma-patterning of NP volume. It would be nice, if we will show %power patterns decencies on intensity for side probe pulse to show beam steering due to symmetry breaking. %\item Fig.4: (a) Dependence on pulse duration is also interesting. We %have to show at which duration the asymmetry factor is saturated. (b) %2D map of asymmetry factor in false colors, where x-axis and y-axis correspond to intensity and NP diameter. %\end{itemize} %Additionally, if you will manage to calculate %evolution of scattering power pattern and show considerable effect of % beam steering, we can try Nanoscale or LPR, because the novelty will % be very high. \section{Modeling details} We focus attention on silicon because this material is promising for the implementation of numerous nonlinear photonic devices. This advantage is based on a broad range of optical nonlinearities, strong two-photon absorption, as well as a possibility of the photo-induced EHP excitation~\cite{leuthold2010nonlinear}. Furthermore, silicon nanoantennas demonstrate a sufficiently high damage threshold due to the high melting temperature ($\approx 1690$~K), whereas its nonlinear optical properties were extensively studied during the last decades~\cite{Van1987, Sokolowski2000, leuthold2010nonlinear}. High silicon melting point typically preserves structures formed from this material up to the EHP densities on the order of the critical value~\cite{Korfiatis2007} $N_{cr} \approx 5\cdot{10}^{21}$~cm$^{-3}$. At the critical density and above, silicon acquires metallic properties ($Re(\epsilon) < 0$) and contributes to the EHP reconfiguration during ultrashort laser irradiation. The process of three-dimensional photo-generation and temporal evolution of EHP in silicon NPs has not been modeled before. Therefore, herein we propose a model considering ultrashort laser interactions with a resonant silicon sphere, where the EHP is generated via one- and two-photon absorption processes. Importantly, we also consider nonlinear feedback of the material by taking into account the intraband light absorption on the generated free carriers. To simplify our model, we neglect free carrier diffusion due to the considered short time scales. In fact, the aim of the present work is to study the EHP dynamics \textit{during} ultra-short (\textit{fs}) laser interaction with the NP. The created electron-hole modifies both laser-particle interaction and, hence, the following particle evolution. However, the plasma then will recombine at picosecond time scale. \subsection{Light propagation} The incident wave propagates in positive direction of $z$ axis, the NP geometric center located at $z=0$ front side corresponds to the volume $z>0$ and back side for $z<0$, as shown in Fig.~\ref{fgr:concept}. Ultra-short laser interaction and light propagation inside the silicon NP are modeled by solving the system of three-dimensional Maxwell's equations written in the following way \begin{align} \begin{cases} \label{Maxwell}$$ \displaystyle{\frac{\partial{\vec{E}}}{\partial{t}}=\frac{\nabla\times\vec{H}}{\epsilon_0\epsilon}-\frac{1}{\epsilon_0\epsilon}(\vec{J}_p+\vec{J}_{Kerr})} \\ \displaystyle{\frac{\partial{\vec{H}}}{\partial{t}}=-\frac{\nabla\times\vec{E}}{\mu_0}}, $$ \end{cases} \end{align} where $\vec{E}$ is the electric field, $\vec{H}$ is the magnetizing field, $\epsilon_0$ is the free space permittivity, $\mu_0$ is the permeability of free space, $\epsilon = n^2 = 3.681^2$ is the permittivity of non-excited silicon at $800$~nm wavelength \cite{Green1995}, $\vec{J}_p$ and $\vec{J}_{Kerr}$ are the nonlinear currents, which include the contribution due to Kerr effect $\vec{J}_{Kerr} = \epsilon_0\epsilon_\infty\chi_3\frac{\partial\left(\left|\vec{E}\right|^2\vec{E}\right)}{\partial{t}}$, where $\chi_3 =4\cdot{10}^{-20}$ m$^2$/V$^2$ for laser wavelength $\lambda = 800$~nm \cite{Bristow2007}, and heating of the conduction band, described by the differential equation derived from the Drude model \begin{equation} \label{Drude} \displaystyle{\frac{\partial{\vec{J_p}}}{\partial{t}} = - \nu_e\vec{J_p} + \frac{e^2N_e(t)}{m_e^*}\vec{E}}, \end{equation} where $e$ is the elementary charge, $m_e^* = 0.18m_e$ is the reduced electron-hole mass \cite{Sokolowski2000}, $N_e(t)$ is the time-dependent free carrier density and $\nu_e = 10^{15}$~s$^{-1}$ is the electron collision frequency \cite{Sokolowski2000}. Silicon NP is surrounded by vacuum, where the light propagation is calculated by Maxwell's equations with $\vec{J} = 0$ and $\epsilon = 1$. The system of Maxwell's equations coupled with electron density equation is solved by the finite-difference numerical method \cite{Rudenko2016}, based on the finite-difference time-domain (FDTD) \cite{Yee1966} and auxiliary-differential methods for dispersive media \cite{Taflove1995}. At the edges of the grid, we apply the absorbing boundary conditions related to convolutional perfectly matched layers to avoid nonphysical reflections \cite{Roden2000}. The initial electric field is introduced as a Gaussian slightly focused beam as follows \begin{align} \begin{aligned} \label{Gaussian} {E_x}(t, r, z) = \frac{w_0}{w(z)}{exp}\left(i\omega{t} - \frac{r^2}{{w(z)}^2} - ikz - ik\frac{r^2}{2R(z)} + i\varsigma(z)\right)\\ \times\;{exp}\left(-\frac{4\ln{2}(t-t_0)^2}{\theta^2}\right), \end{aligned} \end{align} where $\theta \approx 130$~\textit{fs} is the temporal pulse width at the half maximum (FWHM), $t_0$ is the time delay, $w_0 = 3{\mu}m$ is the beam waist, $w(z) = {w_0}\sqrt{1+(\frac{z}{z_R})^2}$ is the Gaussian's beam spot size, $\omega = 2{\pi}c/{\lambda}$ is the angular frequency, $\lambda = 800$~nm is the laser wavelength in air, $c$ is the speed of light, ${z_R} = \frac{\pi{{w_0}^2}n_0}{\lambda}$ is the Rayleigh length, $r = \sqrt{x^2 + y^2}$ is the radial distance from the beam's waist, $R_z = z\left(1+(\frac{z_R}{z})^2\right)$ is the radius of curvature of the wavelength comprising the beam, and $\varsigma(z) = {\arctan}(\frac{z}{z_R}) $ is the Gouy phase shift. \subsection{Material ionization} To account for the material ionization induced by a sufficiently intense laser field inside the particle, we couple Maxwell's equations with the kinetic equation for the electron-hole plasma as described below. The time-dependent conduction-band carrier density evolution is described by the rate equation proposed by van Driel \cite{Van1987}. This equation takes into account such processes as photoionization, avalanche ionization and Auger recombination, and is written as \begin{equation} \label{Dens} \displaystyle{\frac{\partial{N_e}}{\partial t} = \frac{N_a-N_e}{N_a}\left(\frac{\sigma_1I}{\hbar\omega} + \frac{\sigma_2I^2}{2\hbar\omega}\right) + \alpha{I}N _e - \frac{C\cdot{N_e}^3}{C\tau_{rec}N_e^2+1},} \end{equation} where $I=\frac{n}{2}\sqrt{\frac{\epsilon_0}{\mu_0}}\left|\vec{E}\right|^2$ is the intensity, $\sigma_1 = 1.021\cdot{10}^3$~cm$^{-1}$ and $\sigma_2 = 0.1\cdot{10}^{-7}$~cm/W are the one-photon and two-photon interband cross-sections \cite{Choi2002, Bristow2007, Derrien2013}, $N_a = 5\cdot{10}^{22}$~cm$^{-3}$ is the particle saturation density \cite{Derrien2013}, $C = 3.8\cdot{10}^{-31}$~cm$^6$/s is the Auger recombination rate \cite{Van1987}, $\tau_{rec} = 6\cdot{10}^{-12}$~s is the minimum Auger recombination time \cite{Yoffa1980}, and $\alpha = 21.2$~cm$^2$/J is the avalanche ionization coefficient \cite{Pronko1998} at the wavelength $800$~nm in air. As we have noted, free carrier diffusion is neglected during and shortly after the laser excitation \cite{Van1987, Sokolowski2000}. In particular, from the Einstein formula $D = k_B T_e \tau/m^* \approx (1$--$\,2)\cdot{10}^{-3}$ m$^2$/s ($k_B$ is the Boltzmann constant, $T_e$ is the electron temperature, $\tau=1$~\textit{fs} is the collision time, $m^* = 0.18 m_e$ is the effective mass), where $T_e \approx 2*{10}^4$ K for $N_e$ close to $N_{cr}$ \cite{Ramer2014}. It means that during the pulse duration ($\approx 130$~\textit{fs}) the diffusion length will be around 10$\,$--15~nm for $N_e$ close to $N_{cr}$. \begin{figure}[ht!] \centering \includegraphics[width=0.495\textwidth]{mie-fdtd-3} \caption{\label{mie-fdtd} (a) First four Lorentz-Mie coefficients ($a_1$, $a_2$, $b_1$, $b_2$) and (b) factors of asymmetry $G_I$, $G_{I^2}$ according to the Mie theory at fixed wavelength $800$~nm. (c, d) Squared intensity distributions at different radii \textit{R} calculated by the Mie theory and (e, f) EHP distribution for low densities $N_e \approx 10^{20}$~cm$^{-3}$ by Maxwell's equations (\ref{Maxwell}, \ref{Drude}) coupled with EHP density equation (\ref{Dens}). (c-f) Incident light propagates from the left to the right along $Z$ axis, electric field polarization $\vec{E}$ is along $X$ axis.} \end{figure} The changes of the real and imaginary parts of the permittivity associated with the time-dependent free carrier response \cite{Sokolowski2000} can be derived from equations (\ref{Maxwell}, \ref{Drude}) and are written as follows \begin{align} \begin{cases} \label{Index} $$ \displaystyle{Re(\epsilon^\ast) = \epsilon -\frac{{e^2}N_e}{\epsilon_0m_e^*(\omega^2+{\nu_e}^2)}} \\ \displaystyle{Im(\epsilon^\ast) = \frac{{e^2}N_e\nu_e}{\epsilon_0m_e^*\omega(\omega^2+{\nu_e}^2)}.} $$ \end{cases} \end{align} \subsection{Mie calculations} A steady-state interaction of a plain electromagnetic wave with a spherical particle has a well-known analytical solution described by a Mie theory~\cite{Bohren1983}. It is only valid in the absence of nonlinear optical response, thus we can compare it against above-mentioned FDTD-EHP model only for small plasma densities, where we can neglect EHP impact to the refractive index. Non-stationary nature of a~\textit{fs} pulse increases the complexity of the analysis. A detailed discussion on the relation between the Mie theory and FDTD-EHP model will be provided in the next section. We used the Scattnlay program to evaluate calculations of Lorentz-Mie coefficients ($a_i$, $b_i$) and near-field distribution~\cite{Ladutenko2017}. This program is available online at GitHub~\cite{Scattnlay-web} under open source license. \section{Results and discussion} \subsection{Asymmetry factors} \begin{figure*}[p] \centering \includegraphics[width=140mm]{time-evolution-I-no-NP.pdf} \caption{\label{time-evolution} Temporal EHP (a, c, e) and asymmetry factor $G_{N_e}$ (b, d, f) evolution for different Si nanoparticle radii of (a, b) $R = 75$~nm, (c, d) $R = 100$~nm, and (e, f) $R = 115$~nm. Pulse duration $130$~\textit{fs} (FWHM). Wavelength $800$~nm in air. (b, d, f) Different stages of EHP evolution shown in Fig.~\ref{plasma-grid} are indicated. The temporal evolution of Gaussian beam intensity is also shown. Peak laser intensity is fixed to be 10$^{12}$~W/cm$^2$. For better visual representation of time scale at a single optical cycle we put a squared electric field profile in all plots in Fig.~\ref{time-evolution} in gray color as a background image (note linear time scale on the left column and logarithmic scale on the right one).} \end{figure*} \begin{figure*} \centering \includegraphics[width=145mm]{plasma-grid.pdf} \caption{\label{plasma-grid} EHP density snapshots inside Si nanoparticle of radii $R = 75$~nm (a-d), $R = 100$~nm (e-h) and $R = 115$~nm (i-l) taken at different times and conditions of excitation (stages $1-4$: (1) first optical cycle, (2) maximum during few optical cycles, (3) quasi-stationary regime, (4) strongly nonlinear regime). $\Delta{Re(\epsilon)}$ indicates the real part change of the local permittivity defined by Equation (\ref{Index}). Pulse duration 130~\textit{fs} (FWHM). Wavelength 800~nm in air. Peak laser intensity is fixed to be 10$^{12}$~W/cm$^2$. } \end{figure*} %\subsection{Effect of the irradiation intensity on EHP generation} We start with a pure electromagnetic problem without EHP generation. We plot in Fig.~\ref{mie-fdtd}(a) Mie coefficients of a Si NP as a function of its size for a fixed laser wavelength $\lambda = 800$~nm. For the NP sizes under consideration most of contribution to the electromagnetic response originates from electric and magnetic dipole (ED and MD), while for sizes near $R \approx 150$~nm a magnetic quadrupole (MQ) response turns to be the main one. The superposition of multipoles defines the distribution of electric field inside the NP. We introduce $G_I$ factor of asymmetry, corresponding to difference between the volume integral of squared electric field in the front side ($I^{front}$) of the NP to that in the back side ($I^{back}$) normalized to their sum: $G_I = (I^{front}-I^{back})/(I^{front}+I^{back})$, where $I^{front}=\int_{(z>0)}|E|^2d{\mathrm{v}}$ and $I^{back}=\int_{(z<0)} |E|^2d{\mathrm{v}}$ are expressed via amplitude of the electric field $|E|$. The factor $G_{I^2}$ was determined in a similar way by using volume integrals of squared intensity to predict EHP asymmetry due to two-photon absorption. Fig.~\ref{mie-fdtd}(b) shows $G$ factors as a function of the NP size. For the NPs of sizes below the first MD resonance, the intensity is enhanced in the front side as in Fig.~\ref{mie-fdtd}(c) and $G_I > 0$. The behavior changes near the size resonance value, corresponding to $R \approx 105$~nm. In contrast, for larger sizes, the intensity is enhanced in the back side of the NP as demonstrated in Fig.~\ref{mie-fdtd}(d). In fact, rather similar EHP distributions can be obtained by applying Maxwell's equations coupled with the rate equation for a relatively weak excitation with EHP concentration of $N_e \approx 10^{20}$~cm$^{-3}$, see Fig.~\ref{mie-fdtd}(e,f). The optical properties do not change considerably due to the excitation according to (\ref{Index}). Therefore, the excitation processes follow the intensity distribution. However, such coincidence was achieved under quasi-stationary conditions, after the electric field made enough oscillations inside the Si NP. Further on, we present transient analysis, which reveals much more details. To achieve a quantitative description for evolution of the EHP distribution during the \textit{fs} pulse, we introduced another asymmetry factor $G_{N_e} = (N_e^{front}-N_e^{back})/(N_e^{front}+N_e^{back})$ indicating the relationship between the average EHP densities in the front and in the back halves of the NP, defined as $N_e^{front}=\frac{2}{V}\int_{(z>0)} {N_e}d{\mathrm{v}}$ and $N_e^{back}=\frac{2}{V}\int_{(z<0)} {N_e}d{\mathrm{v}}$, where $V = \frac{4}{3}\pi{R}^3$ is the volume of the nanosphere. In this way, $G_{N_e} = 0$ corresponds to the symmetrical case and the assumption of the NP homogeneous EHP distribution can be made to investigate the optical response of the excited Si NP. When $G_{N_e}$ significantly differs from $0$, this assumption, however, could not be justified. In what follows, we discuss the results of the numerical modeling of the temporal evolution of EHP densities and the asymmetry factors $G_{N_e}$ for different sizes of Si NP, as shown in Fig.~\ref{time-evolution}). Typical change in permittivity corresponding to each stage is shown in Fig.~\ref{plasma-grid}. It reveals the EHP evolution stages during interaction of femtosecond laser pulse with the Si NPs. \subsection{Stages of transient Si nanoparticle photoexcitation} To describe all the stages of non-linear light interaction with Si NP, we present the calculation results obtained by using Maxwell's equations coupled with electron kinetics equations for different radii for resonant and non-resonant conditions at peak intensity 10$^{12}$~W/cm$^2$ and $\lambda = 800~nm$. In this case, the geometry of the EHP distribution can strongly deviate from the intensity distribution given by the Mie theory. Two main reasons cause the deviation: (i) non-stationarity of interaction between electromagnetic pulse and NP (ii) nonlinear effects, taking place due to transient optical changes in Si. The non-stationary intensity deposition during \textit{fs} pulse results in different time delays for exciting electric and magnetic resonances inside Si NP because of different quality factors $Q$ of the resonances. In particular, MD resonance (\textit{b1}) has $Q \approx 8$, whereas electric one (\textit{a1}) has $Q \approx 4$. The larger particle supporting MQ resonance (\textit{b2}) demonstrates $ Q \approx 40$. As soon as the electromagnetic wave period at $\lambda = 800$~nm is $\approx 2.7$~\textit{fs}, one needs about 10~\textit{fs} to pump the ED, 20~\textit{fs} for the MD, and about 100~\textit{fs} for the MQ. According to these considerations, after few optical cycles taking place on a 10~\textit{fs} scale it results in the excitation of the low-\textit{Q} ED resonance, which dominates MD and MQ independently on the exact size of NPs. Moreover, during the first optical cycle there is no multipole modes structure inside NP, which results in a very similar field distribution for all sizes of NP under consideration as shown in Fig.~\ref{plasma-grid}(a,e,i). We address to this phenomena as \textit{'Stage~1'}. This stage demonstrates the initial penetration of electromagnetic field into the NP during the first optical cycle. Resulting factors $G_{N_e}$ exhibit sharp increase at this stage (Fig.~\ref{time-evolution}(b,d,f), yielding strong and ultrafast symmetry breaking. \textit{'Stage~2'} corresponds to further electric field oscillations ($t \approx 5$--$15$~\textit{fs}) leading to the formation of ED E-field pattern in the center of the Si NP as shown in Fig.~\ref{plasma-grid}(f,j). We stress the nonstationary nature of E-field pattern at this stage, whereas there is a simultaneous growth of the incident pulse amplitude. This leads to a superposition of ED near-field pattern with that from Stage 1, resulting in EHP concentration in the front side of the Si NP. This effect dominates for the smallest NP with $R=75$~nm in Fig.~\ref{plasma-grid}(b), where ED mode is tuned far away from the resonance (see Fig.~\ref{mie-fdtd}(c) for field suppression inside NP predicted by the Mie theory). At this stage, the density of EHP ($N_e < 10^{20}$~cm$^2$) is still not high enough to significantly affect the optical properties of the NP (Figs.~\ref{time-evolution}(a,c,e)). When the number of optical cycles is large enough ($t>20$~\textit{fs}) both ED and MD modes can be exited to the level necessary to achieve the stationary intensity pattern corresponding to the Mie-based intensity distribution at \textit{'Stage~3'} (see Fig.~\ref{plasma-grid}(c,g,k)). The EHP density for the most volume of NP is still relatively small to affect the EHP evolution, but is already high enough to change the local optical properties, i.e. real part of permittivity. Below the MD resonance ($R = 75$~nm), the EHP is mostly localized in the front side of the NP as shown in Fig.~\ref{plasma-grid}(c). The highest quasi-stationary asymmetry factor $G_{N_e} \approx 0.5$--$0.6$ is achieved in this case (Fig~\ref{time-evolution}(b)). At the MD resonance conditions ($R = 100$~nm), the EHP distribution has a toroidal shape and is much closer to the homogeneous distribution. In contrast, above the MD resonant size for $R = 115$~nm the $G_{N_e} < 0$ due to the fact that EHP is dominantly localized in the back side of the NP. Diue to a quasi-stationary pumping during Stage~3, it is superposed with Stage~1 field pattern, resulting in an additional EHP localized in the front side. This can be seen when comparing result from the Mie theory in Fig.~\ref{mie-fdtd}(d) with that of full 3D simulation in Fig.~\ref{mie-fdtd}(f). Note that pumping of NP significantly changes during a single optical cycle, this leads to a large variation of asymmetry factor $G_{N_e}$ at the first stage. This variation steadily decreases as it goes to Stage~3, as shown in Fig.~\ref{time-evolution}(b,d,f). To explain this effect, we consider the time evolution of average EHP densities $N_e$ in the front and back halves of the NP presented in Fig.~\ref{time-evolution}(a,c,e). As soon as the recombination and diffusion processes are negligible at \textit{fs} time scale, both $N_e^{front}$ and $N_e^{back}$ curves experience monotonous behavior with small pumping steps synced to the incident pulse. The front and the back halves of NP are separated in space, which obviously leads to the presence of time delay between pumping steps in each curve caused by the same optical cycle of the incident wave. This delay causes a large asymmetry factor during the first stage. However, as soon as average EHP density increases, the relative contribution of these pumping steps to the resulting asymmetry becomes smaller. This way variations of asymmetry $G_{N_e}$ synced with the period of incident light decrease. Higher excitation conditions are followed by larger values of electric field amplitude, which leads to the appearance of high EHP densities causing a significant change in the optical properties of silicon according to the equations (\ref{Index}). From the Mie theory, the initial (at the end of Stage~3) spatial pattern of the optical properties is non-homogeneous. When non-homogeneity of the optical properties becomes strong enough it leads to the reconfiguration of the E-field inside NP, which in turn strongly affects further reconfiguration of the optical properties. We refer to these strong nonlinear phenomena as \textit{'Stage~4'}. In general, the reconfiguration of the electric field is unavoidable as far as the result from the Mie theory comes with the assumption of homogeneous optical properties in a spherical NP. Thus, the evolution of EHP density during Stage~4 depends on the result of multipole modes superposition at the end of Stage~3 and is quite different as we change the size of NP. For $R=75$~nm and $R=100$~nm, we observe a front side asymmetry before Stage~4, however, its origin is quite different. The $R=75$~nm NP is out of resonance, moreover, Mie field pattern and the one, which comes from Stage~1 are quite similar. As soon as EHP density becomes high enough to change optical properties, the NP is still out of resonance, however, the presence of EHP increases absorption in agreement with (\ref{Index}), decreases Q-factor, and destroys optical modes. % This effect effectively leads to a partial screening, and it % becomes harder for the incident wave to penetrate deeper into % EHP. Finally, this finishes spilling the NP`s volume with plasma % reducing the asymmetry, see Fig.~\ref{plasma-grid}(d). For $R=100$~nm, the evolution during the final stage goes in a similar way, with a notable exception regarding MD resonance. As soon as presence of EHP increases the absorption, it suppresses the MD resonance with symmetric filed pattern, thus, the asymmetry factor can be increased. This result was observed in Fig.~\ref{time-evolution}(d) with a local maximum near 100~\textit{fs} mark. The last NP with $R=115$~nm shows the most complex behavior during Stage~4. The superposition of Mie-like E-field pattern with that from Stage~1 results in the presence of two EHP spatial maxima, back and front shifted. They serve as starting seeds for the EHP formation, and an interplay between them forms a complex behavior of the asymmetry factor curve. Namely, the sign is changed from negative to positive and back during the last stage. This numerical result can hardly be explained in a simple qualitative manner, it is too complex to account all near-field interaction of incident light with two EHP regions inside a single NP. It is interesting to note, however, that in a similar way as it was for $R=100$~nm the increased absorption should destroy ED and MD resonances, which are responsible for the back-shifted EHP. As soon as this EHP region is quite visible on the last snapshot in Fig.~\ref{plasma-grid}(l), this means that EHP seeds are self-supporting. %A bookmark by Kostya As the EHP acquires quasi-metallic properties at stronger excitation $N_e > 5\cdot{10}^{21}$~cm$^{-3}$, the EHP distribution evolves inside NPs because of the photoionization and avalanche ionization induced transient optical response and the effect of newly formed EHP. This way, the distribution becomes more homogeneous and the effect is likely to be enhanced by electron diffusion inside Si NPs. It is worth noting that it is possible to achieve the formation of deeply subwavelength EHP regions due to high field localization. The smallest EHP localization and the larger asymmetry factor are achieved below the MD resonant conditions for $R < 100$~nm. Thus, the EHP distribution in Fig.~\ref{plasma-grid}(c) is optimal for symmetry breaking in Si NP, as it results in the larger asymmetry factor $G_{N_e}$ and higher electron densities $N_e$. We stress here that such regime could be still safe for NP due to the very small volume where such high EHP density is formed. % \subsection{Effects of NP size and scattering efficiency % factor on scattering directions} % \begin{figure}[ht] \centering % \includegraphics[width=90mm]{time-evolution.png} % \caption{\label{time-evolution} a) Scattering efficiency factor $Q_{sca}$ % dependence on the radius $R$ of non-excited silicon NP % calculated by Mie theory; b) Parameter of forward/backward scattering % dependence on the radius $R$ calculated by Mie theory for non-excited % silicon NP c) Optimization parameter $K$ dependence on the % average electron density $n_e^{front}$ in the front half of the % NP for indicated radii (1-7).} % \end{figure} % We have discussed the EHP kinetics for a silicon NP of a % fixed radius $R \approx 105$~nm. In what follows, we investigate the % influence of the NP size on the EHP patterns and temporal % evolution during ultrashort laser irradiation. A brief analysis of % the initial intensity distribution inside the NP given by % the classical Mie theory for homogeneous spherical particles % \cite{Mie1908} can be useful in this case. Fig. \ref{time-evolution}(a, b) % shows the scattering efficiency and the asymmetry parameter for % forward/backward scattering for non-excited silicon NPs of % different radii calculated by Mie theory \cite{Mie1908}. Scattering % efficiency dependence gives us the value of resonant sizes of % NPs, where the initial electric fields are significantly % enhanced and, therefore, we can expect that the following conditions % will result in a stronger electron density gradients. Additionally, % in the case of maximum forward or backward scattering, the initial % intensity distribution has the maximum of asymmetry. One can note, % that for $R \approx 100$~nm and $R \approx 150$~nm both criteria are % fulfilled: the intensity is enhanced $5-10$ times due to % near-resonance conditions and its distribution has a strong % asymmetry. % In what follows, we present the calculation results obtained by % using Maxwell's equations coupled with electron kinetics for % different extremum radii for resonant and non-resonant % conditions. One can note, that the maximum asymmetry factor of EHP % $G$ does not guarantee the optimal asymmetry of intensity % distribution, as the size of generated plasma and the value of the % electron density equally contribute to the change of the modified % NP optical response. For example, it is easier to localize % high electron densities inside smaller NPs, however, due % to the negligible size of the generated EHP with respect to laser % wavelength in media, the intensity distribution around the % NP will not change considerably. Therefore, we propose to % introduce the optimization factor % $K = \frac{n_e(G-1)R^2}{n_{cr}G{R_0^2}}$, where $R_0 = 100$~nm, % $n_{cr} = 5\cdot{10}^{21}~cm^{-3}$, and $G$ is asymmetry of EHP, % defined previously. The calculation results for different radii of % silicon NPs and electron densities are presented in % Fig. \ref{time-evolution}(c). One can see, that the maximum value are achieved % for the NPs, that satisfy both initial maximum forward % scattering and not far from the first resonant condition. For larger % NPs, lower values of EHP asymmetry factor are obtained, as % the electron density evolves not only from the intensity patterns in % the front side of the NP but also in the back side. %TODO: %Need to discuss agreement/differences between Mie and FDTD+rate. Anton ?! % To demonstrate the effect of symmetry breaking, we calculate the % intensity distribution around the NP for double-pulse % experiment. The first pulse of larger pulse energy and polarization % along $Ox$ generates asymmetric EHP inside silicon NP, % whereas the second pulse of lower pulse energy and polarization $Oz$ % interacts with EHP after the first pulse is gone. The minimum % relaxation time of high electron density in silicon is % $\tau_{rec} = 6\cdot{10}^{-12}$~s \cite{Yoffa1980}, therefore, the % electron density will not have time to decrease significantly for % subpicosecond pulse separations. In our simulations, we use % $\delta{t} = 200\:f\!s$ pulse separation. The intensity % distributions near the silicon NP of $R = 95$~nm, % corresponding to maxima value of $K$ optimization factor, without % plasma and with generated plasma are shown in Fig. \ref{fig4}. The % intensity distribution is strongly asymmetric in the case of EHP % presence. One can note, that the excited NP is out of % quasi-resonant condition and the intensity enhancements in % Fig. \ref{fig4}(c) are weaker than in Fig. \ref{fig4}(b). Therefore, % the generated nanoplasma acts like a quasi-metallic nonconcentric % nanoshell inside the NP, providing a symmetry reduction % \cite{Wang2006}. % \begin{figure}[ht] \centering % \includegraphics[width=90mm]{fig4.png} % \caption{\label{fig4} a) Electron plasma distribution inside Si % NP $R \approx 95$~nm 50~\textit{fs} after the pulse peak; % (Kostya: Is it scattering field intensity snapshot $XXX\:f\!s$ after % the second pulse maxima passed the particle?) Intensity % distributions around and inside the NP b) without plasma, % c) with electron plasma inside.} % \end{figure} %\begin{figure} %\centering % \includegraphics[height=0.7\linewidth]{Si-flow-R140-XYZ-Eabs} % \caption{EHP distributions for nonres., MD, ED, and MQ NPs % at moderate photoexcitation. The aim is to show different possible % EHP patterns and how strong could be symmetry breaking. % \label{fgr:example} %\end{figure} %\subsection{Asymmetry analysis: effects of pulse duration, intensity % and size} It is important to optimize asymmetry by varying pulse % duration, intensity and size. \section{Conclusions} We have rigorously modeled and studied ultra-fast and intense light interaction with a single silicon nanoparticle of various sizes for the first time to our best knowledge. As a result of the presented self-consistent nonlinear calculations, we have obtained spatio-temporal EHP evolution inside the NPs and investigated the asymmetry of the EHP distributions. We have revealed EHP strong asymmetric distribution during the first optical cycle for different sizes. The highest average EHP asymmetry has been observed for NPs of smaller sizes below the first magnetic dipole resonance, when EHP is concentrated in the front side mostly during the laser pulse absorption. Essentially different EHP evolution and lower asymmetry has been achieved for larger NPs due to the intensity enhancement in the back side of the NP. The EHP densities above the critical value have been shown to lead to homogenization of the EHP distribution. The observed plasma-induced breaking symmetry can be useful for creation of nonsymmetrical nanophotonic designes, e.g. for beam steering or enhanced second harmonics generation. Also, the asymmetric EHP opens a wide range of applications in NP nanomashining at deeply subwavelength scale. \section{Acknowledgments} A. R. and T. E. I. gratefully acknowledge the CINES of CNRS for computer support. S. V. M. is thankful to ITMO Fellowship Program. This work was partially supported by Russian Foundation for Basic Researches (grants 17-03-00621, 17-02-00538, 16-29-05317). %%%END OF MAIN TEXT%%% %The \balance command can be used to balance the columns on the final %page if desired. It should be placed anywhere within the first column %of the last page. %\balance %If notes are included in your references you can change the title % from 'References' to 'Notes and references' using the following % command: % \renewcommand\refname{Notes and references} %%%REFERENCES%%% \bibliography{References} %You need to replace "rsc" on this line %with the name of your .bib file \bibliographystyle{rsc} %the RSC's .bst file \end{document}