## 2 Methods
The MNP is described as a continuum of energy states, while the asymmetric double SQD molecule is depicted as a three-level \(\Lambda\)-type system with a ground state \(|0\rangle\), a direct exciton state \(|1\rangle\) (electron-hole pair in the same quantum dot), and an indirect exciton state \(|2\rangle\) (due to the transferring of the excited electron to the excited state of the second quantum dot). The application of a threshold voltage to the double SQD induces strong tunneling between the upper levels, while, at the same time, the ground energy levels are effectively decoupled. where  
\[\Omega_{i}=\frac{E_{i}\mu}{2\hbar e_{eff}s}\bigg{(}1+\frac{s_{a}\gamma_{1}a^{3}}{R^ {3}}\bigg{)}, \tag{3}\]  
\[G_{i}=\sum_{n=1}^{N}\frac{1}{4\pi\epsilon_{env}}\frac{s_{n}\gamma_{n}(\omega_{i} )a^{2n+1}\mu^{2}}{\hbar e_{eff}^{2}S^{2n+4}}, \tag{4}\]  
with \(\gamma_{n}(\omega_{i})=[-\varepsilon_{env}+\varepsilon_{m}(\omega_{i})]/[(n+ 1)\varepsilon_{env}/n+\varepsilon_{m}(\omega_{i})]\) and \(\varepsilon_{eff}=(\varepsilon_{S}+2\varepsilon_{env})/3\varepsilon_{env}\). In the above, we have used the quasistatic approximation, where retardation effects are omitted, following several previous works [2; 3; 4; 5; 6; 7; 8; 9; 10; 11; 12; 13; 14; 15; 16; 17; 18; 19; 20; 21; 22; 23; 24; 25; 26; 27; 28; 29; 30; 31; 32; 33; 34; 35; 36; 37; 38; 39; 40; 41; 42; 43; 44; 45; 46; 47; 48; 49; 53; 54; 55; 56]. The first term in the expression of Equation (3) is ascribed to the direct coupling between the applied electromagnetic field and the SQDs. We expect that the MNP interacts with the applied field, as well as the field emitted by the SQDs, and as a result, it emits in its surrounding space. The coupling between the SQDs and the first (second) component of this emitted field is associated with the second term of Equation (3) (the \(G_{i}\) parameter). The parameters \(s_{a}\) and \(s_{n}\), in the case of an incident field polarized along the direction of the interparticle axis, are respectively equal to 2 and \((n+1)^{2}\), while, for an incident field that is polarized perpendicularly to the direction of the interparticle axis, we take \(s_{a}=-1\) and \(s_{n}=n(n+1)/2\). The \(G_{i}\) parameter is calculated including the multipole effects for the MNP [55]. Based on the expression of the Hamiltonian that describes the system, in the rotating wave approximation, a set of nine density-matrix differential equations is derived. The number of these equations is eliminated to five, under the conditions \(\sum_{0}^{2}\rho_{ii}=1\) and \(\rho_{mn}=\rho_{nm}^{*}\), with \(n\
eq m\):  
\[ \begin{array}{rl}\dot{\rho}_{01}(t)=&i(\Delta+i\gamma_{01})\rho_{01}(t)+iT_ {e}\rho_{02}(t)\\ &-i\big{[}\Omega_{a}^{*}+G_{a}^{*}\rho_{01}(t)+\Omega_{b}^{*}e^{i\beta t}+G_{b }^{*}\rho_{01}(t)\big{]}\big{[}\rho_{00}(t)-\rho_{11}(t)\big{]},\end{array}  \tag{5}\]  
\[ \begin{array}{rl}\dot{\rho}_{02}(t)=&-i(-\Delta+\omega_{12}-i\gamma_{02}) \rho_{02}(t)+iT_{e}\rho_{01}(t)\\ &+i\big{[}\Omega_{a}^{*}+G_{a}^{*}\rho_{01}(t)+\Omega_{b}^{*}e^{i\beta t}+G_{b }^{*}\rho_{01}(t)\big{]}\rho_{12}(t),\end{array}  \tag{6}\]  
\[ \begin{array}{rl}\dot{\rho}_{12}(t)=&i(-\omega_{12}+i\gamma_{12})\rho_{12}(t )-iT_{e}\big{[}1-\rho_{00}(t)-2\rho_{11}(t)\big{]}\\ &+i\big{[}\Omega_{a}+G_{a}\rho_{10}(t)+\Omega_{b}e^{-i\beta t}+G_{b}\rho_{10}( t)\big{]}\rho_{02}(t),\end{array}  \tag{7}\]  
\[ \begin{array}{rl}\dot{\rho}_{00}(t)=&-i\big{[}\Omega_{a}+G_{a}\rho_{10}(t)+ \Omega_{b}e^{-i\beta t}+G_{b}\rho_{10}(t)\big{]}\rho_{01}(t)+i[\Omega_{a}^{*}+ G_{a}^{*}\rho_{01}(t)\\ &+\Omega_{b}^{*}e^{i\beta t}+G_{b}^{*}\rho_{01}(t)\big{]}\rho_{10}(t)-\Gamma_{ 20}[\rho_{00}(t)+\rho_{11}(t)-1]+\Gamma_{10}\rho_{11}(t),\end{array}  \tag{8}\]  
\[ \begin{array}{rl}\dot{\rho}_{11}(t)=&i\Big{[}\Omega_{a}+G_{a}\rho_{10}(t)+ \Omega_{b}e^{-i\beta t}+G_{b}\rho_{10}(t)\Big{]}\rho_{01}(t)-i[\Omega_{a}^{*}+ G_{a}^{*}\rho_{01}(t)\\ &+\Omega_{b}^{*}e^{i\beta t}+G_{b}^{*}\rho_{01}(t)\Big{]}\rho_{10}(t)+iT_{e}[ \rho_{12}(t)-\rho_{21}(t)]-(\Gamma_{10}+\Gamma_{12})\rho_{11}(t),\end{array}