Mode Coupling Properties of the Plasmonic Dimers Composed of Graphene Nanodisks

The electromagnetic properties of the plasmonic dimer composed of coupled graphene nanodisks are numerically investigated in this paper. The results demonstrate that the degeneracy of the plasmonic modes of the dimer is lifted when the coupling is introduced. The evolution of the plasmonic mode, with the variation of inter-disk distance and the chemical potential of one of the nanodisks, is studied. The proposed structure might find broad areas of application including light–matter interaction, optical switching, directional emission of the plasmonic emitter, etc.


Introduction
Similar to electronic atom and molecule systems, photonic molecules (PhMs) and plasmonic molecules (PMs) are composed of interacting nanostructures in which the specific collective electromagnetic (EM) behaviors are demonstrated [1,2]. These behaviors include Fano resonance [3][4][5], slow light [6,7], EM-induced transparency [8][9][10] etc., which result from the coupling between the plasmons and the external disturbance. In the past decade, PhMs and PMs have found significant applications in broad areas including fundamental sciences, the processing of information and biochemics due to their excellent properties of high quality (Q) factors, flexible structure design, and high sensitivity to the variation of the environments [11]. Among the reported works, Qiang Li et al. [12] proposed and fabricated a hybrid photonic-plasmonic molecule composed of two coupled identical Ag-Al 2 O 3 -Si sandwich disks. This hybrid photonic-plasmonic molecule has sub-wavelength properties and provides a novel platform for the exploration of the light-matter interaction. Sang Wook Kim et al. [13] reported the energy-level evolution of coupled non-identical microdisks (CNMs), which were characterized by their avoided crossing which is intrinsically associated with the exceptional point or the non-Hermitian degeneracy. Svetlana V. Boriskina [14] analyzed the wavelength of the Q factors of the whispering-gallery (WG) modes in the PhM micro-disks as a function of the parameter of the PhM numerically, and demonstrated the crossing and avoid-crossing characteristics. However, all the above-mentioned micro or nanodisk based cavities were made by a conventional novel metal such as Ag and Au, or by dielectric materials such as SiO 2 . Thus, the huge ohmic loss would degrade the properties of the devices, for example, Q factors. Furthermore, frequency tunability is impossible as long as the geometric structures of the devices are determined. All the intrinsic shortcomings of the conventional plasmonic materials become an obstacle to practical application.
Recently, graphene, which is composed of a single layer of carbon atoms arranged in a two-dimensional (2D) honeycomb lattice, has attracted great attention in the fields of electronics and photonics all over the world due to its excellent material properties [15][16][17]. More specifically, in the field of the plasmonics, graphene offers the advantages of high confinement capability, low damping loss and tunability of working frequency. Consequently, graphene is a reliable material to substitute conventional metal materials in the field of plasmonics.
Zhu X. et al. investigated plasmon-phonon coupling in large-area graphene dot and antidot arrays [18]. Rosolen G. and Maes B. reported asymmetric and connected graphene dimers for a tunable plasmonic response [19]. Hadad Y. et al. showed that an electrically small chemical potential interruption in a graphene layer may possess a quasi-static resonance that launches a mixed-type semilocalized plasmon wave [20]. They also designed a system that is extremely sensitive to the quantized spectrum of graphene Landau levels and predicted up to 90 • of Faraday rotation with just a single sheet of graphene [21]. We have proposed and numerically explored the plasmonic dimers composed of coupled graphene nanodisks of different sizes, which are surrounded by another chemical potential on the same piece of infinite graphene monolayer [22]. The split modes with non-degenerate wavelengths appear in the proposed structure. Further, the coupling between fundamental and first order WG modes has been pointed out. In the current work, we further analyze the field properties of the collective mode including the resonant frequency; the Q factor in a graphene plasmonic dimer with different chemical potentials in each nanodisk; and we discuss the physical principles behind the EM behaviors of the dimers. The proposed plasmonic molecules could find broad areas of application in fields including tunable filters [23], optical sensors [24] and optical switchers [25].

Models and Simulation Method
The 2D schematic of the plasmonic dimer molecule composed of graphene nanodisks is illustrated in Figure 1, where R is the radius of the nanodisks; g is the gap between the two nanodisks; the chemical potential of each nanodisk is µ c2 and µ c3 , respectively. The chemical potential of the infinite area is µ c1 . Perfectly matched layer (PML) is employed as the boundary condition in order to avoid the reflected wave at the boundaries. In this article, the plasmonic properties of the dimers are investigated by utilizing the commercially available software COMSOL Multi-physics, RF module, eigen frequency solver. The graphene monolayer is modeled as an ultrathin film which can be characterized by the effective index defined as n eff = β/k 0 , where k 0 = 2π/λ is the wavenumber in free space [26]. The propagation constant of surface plasmon polaritons (SPPs) supported by monolayer graphene can be written as [27,28] where η 0 (377 Ω) is the intrinsic impedance of the air, which can be expressed by (µ 0 /ε 0 ) 1/2 where µ 0 and ε 0 represent the permeability and permittivity of vacuum, respectively. The surface conductivity of graphene σ g , composed of the interband electron transitions σ inter and the intraband electron-photon scattering σ intra , is given by the Kubo formula [17,[27][28][29], where µ c is the chemical potential of graphene; ω is the angular frequency of the plasmon; is the reduced Planck's constant; e is the electron charge; k B is the Boltzmann's constant; T is the temperature; and τ is the electron momentum relaxation time. Specifically, the chemical potential of graphene is tunable by either chemical doping or external electrical gating locally [27,30]. Pioneering work showed that a chemical potential as high as 2 eV has been obtained experimentally [31]. Also, a relaxation time as high as 3ps has been achieved experimentally [32]. Thus, in this paper, a relaxation time less than 1ps and a chemical potential less than 1 eV are chosen, which are conservative enough to ensure the reliability of our calculation.
Appl. Sci. 2017, 7, 359 3 of 10 where μc is the chemical potential of graphene; ω is the angular frequency of the plasmon; ћ is the reduced Planck's constant; e is the electron charge; kB is the Boltzmann's constant; T is the temperature; and τ is the electron momentum relaxation time. Specifically, the chemical potential of graphene is tunable by either chemical doping or external electrical gating locally [27,30]. Pioneering work showed that a chemical potential as high as 2 eV has been obtained experimentally [31]. Also, a relaxation time as high as 3ps has been achieved experimentally [32]. Thus, in this paper, a relaxation time less than 1ps and a chemical potential less than 1 eV are chosen, which are conservative enough to ensure the reliability of our calculation.

Results and Discussion
When two nanodisks approach each other, the EM fields of the plasmon interact, and the collective modes appear. First, we study the collective resonant EM mode field distribution as a function of the gap between the two nanodisks, where μc2 = μc3 = 0.5 eV, μc1 = 0.1 eV and R = 50 nm. The z component of the typical electric field distribution of TM2,1 modes of the dimer with a gap of 0 nm and 15 nm are illustrated in Figure 2a,b, respectively. According to the symmetry with respect to the X and Y directions, the four categories of the graphene plasmonic dimer modes are defined as [|+1>x,|+1>y], which is symmetric with respect to both the X and Y directions; [|+1>x,|−1>y], which is symmetric with respect to the X direction and anti-symmetric with respect to the Y direction; [|−1>x,|+1>y], which is anti-symmetric with respect to the X direction and symmetric with respect to the Y direction; and [|−1>x,|−1>y], which is anti-symmetric with respect to both the X and Y direction. More specifically, we follow Ref. [12] and regard the modes which are anti-symmetric with respect to the Y direction as the anti-bonding modes. Normally, the anti-bonding modes are located at the higher energy states. On the other hand, the modes that are symmetric with respect to the Y direction are regarded as the bonding modes, and are located at the lower energy states [12]. The frequency and the quality (Q) factor of the collective modes, as a function of the gap g, are plotted in Figure 3a,b, respectively. It is shown in Figure 3a that the frequencies of the anti-bonding modes decrease with the increasing of the gap g, while the counterparts of the bonding modes demonstrate the converse trend. The coupling strength between the WG modes of the nanodisks is indicated by the frequency splitting of the bonding or anti-bonding modes. When the gap g increases, the coupling strength between the plasmonic modes of the individual nanodisks degrades. Thus, the frequency splitting of the anti-bonding and the bonding modes decreases accordingly [12]. When g further increases to more than 50 nm, all the four mode categories combine together, the result of which is identical to the resonant mode of a single isolated nanodisk. This means that there is no interaction between the two nanodisks of the dimer.

Results and Discussion
When two nanodisks approach each other, the EM fields of the plasmon interact, and the collective modes appear. First, we study the collective resonant EM mode field distribution as a function of the gap between the two nanodisks, where µ c2 = µ c3 = 0.5 eV, µ c1 = 0.1 eV and R = 50 nm. The z component of the typical electric field distribution of TM 2,1 modes of the dimer with a gap of 0 nm and 15 nm are illustrated in Figure 2a,b, respectively. According to the symmetry with respect to the X and Y directions, the four categories of the graphene plasmonic dimer modes are defined as [|+1> x ,|+1> y ], which is symmetric with respect to both the X and Y directions; [|+1> x ,|−1> y ], which is symmetric with respect to the X direction and anti-symmetric with respect to the Y direction; [|−1> x ,|+1> y ], which is anti-symmetric with respect to the X direction and symmetric with respect to the Y direction; and [|−1> x ,|−1> y ], which is anti-symmetric with respect to both the X and Y direction. More specifically, we follow Ref. [12] and regard the modes which are anti-symmetric with respect to the Y direction as the anti-bonding modes. Normally, the anti-bonding modes are located at the higher energy states. On the other hand, the modes that are symmetric with respect to the Y direction are regarded as the bonding modes, and are located at the lower energy states [12]. The frequency and the quality (Q) factor of the collective modes, as a function of the gap g, are plotted in Figure 3a,b, respectively. It is shown in Figure 3a that the frequencies of the anti-bonding modes decrease with the increasing of the gap g, while the counterparts of the bonding modes demonstrate the converse trend. The coupling strength between the WG modes of the nanodisks is indicated by the frequency splitting of the bonding or anti-bonding modes. When the gap g increases, the coupling strength between the plasmonic modes of the individual nanodisks degrades. Thus, the frequency splitting of the anti-bonding and the bonding modes decreases accordingly [12]. When g further increases to more than 50 nm, all the four mode categories combine together, the result of which is identical to the resonant mode of a single isolated nanodisk. This means that there is no interaction between the two nanodisks of the dimer.  It is also shown in Figure 3b that the Q factor of the collective plasmonic mode of the graphene dimer demonstrates an identical trend to that of the resonant frequency as a function of the gap g.
To understand this behavior, we plot the real part and imaginary part of the propagation constant β as a function of frequency in Figure 4, where the chemical potential is 0.5 eV. In general, the Q factor is given by [33,34], where rad Q relates to the radiation loss and abs Q relates to the absorption loss. The radiation loss is determined by the radius and material of the cavity, and abs Q is calculated by where Re(β) and |Im(β)| are the real part and the absolute value of the imaginary part of the propagation constant β. For a graphene sheet with a chemical potential of 0.5 eV, when the frequency of the plasmon increases from 80 THz to 83.5 THz, Re(β) and |Im(β)| increase together. However, the rising rate of Re(β) is even much faster than 2|Im(β)|. Consequently, the Q factor increases as the increasing of frequency, and finally results in a similar frequency trend as a function of given g, which are shown in Figure 3a,b. Next, we investigate the frequency splitting ΔF of the corresponding anti-bonding and bonding modes as a function of the azimuthal mode number m of the plasmonic WG mode TMm,1 of the graphene dimers. Figure 5 shows the frequency splitting ΔF of bonding and anti-bonding of WG modes as a function of m, where the distance g is 5 nm. When the azimuthal number increases from 2 to 8, the frequency ΔF of the bonding mode decreases from 0.567 to 0.002 THz, while the counterpart of the anti-bonding mode decreases from 0.275 to 0.0002 THz. When the azimuthal number increases, the wavelength reduces. However, the EM field confinement capacity of the graphene dimer becomes stronger as the resonant wavelength decreases. So the coupling between the plasmonic modes of the nanodisks becomes weaker, which further results in the small frequency splitting for both the bonding and anti-bonding mode, as shown in Figure 5.  It is also shown in Figure 3b that the Q factor of the collective plasmonic mode of the graphene dimer demonstrates an identical trend to that of the resonant frequency as a function of the gap g.
To understand this behavior, we plot the real part and imaginary part of the propagation constant β as a function of frequency in Figure 4, where the chemical potential is 0.5 eV. In general, the Q factor is given by [33,34], where rad Q relates to the radiation loss and abs Q relates to the absorption loss. The radiation loss is determined by the radius and material of the cavity, and abs Q is calculated by where Re(β) and |Im(β)| are the real part and the absolute value of the imaginary part of the propagation constant β. For a graphene sheet with a chemical potential of 0.5 eV, when the frequency of the plasmon increases from 80 THz to 83.5 THz, Re(β) and |Im(β)| increase together. However, the rising rate of Re(β) is even much faster than 2|Im(β)|. Consequently, the Q factor increases as the increasing of frequency, and finally results in a similar frequency trend as a function of given g, which are shown in Figure 3a,b. Next, we investigate the frequency splitting ΔF of the corresponding anti-bonding and bonding modes as a function of the azimuthal mode number m of the plasmonic WG mode TMm,1 of the graphene dimers. Figure 5 shows the frequency splitting ΔF of bonding and anti-bonding of WG modes as a function of m, where the distance g is 5 nm. When the azimuthal number increases from 2 to 8, the frequency ΔF of the bonding mode decreases from 0.567 to 0.002 THz, while the counterpart of the anti-bonding mode decreases from 0.275 to 0.0002 THz. When the azimuthal number increases, the wavelength reduces. However, the EM field confinement capacity of the graphene dimer becomes stronger as the resonant wavelength decreases. So the coupling between the plasmonic modes of the nanodisks becomes weaker, which further results in the small frequency splitting for both the bonding and anti-bonding mode, as shown in Figure 5. It is also shown in Figure 3b that the Q factor of the collective plasmonic mode of the graphene dimer demonstrates an identical trend to that of the resonant frequency as a function of the gap g.
To understand this behavior, we plot the real part and imaginary part of the propagation constant β as a function of frequency in Figure 4, where the chemical potential is 0.5 eV. In general, the Q factor is given by Q −1 total = Q −1 rad + Q −1 abs [33,34], where Q rad relates to the radiation loss and Q abs relates to the absorption loss. The radiation loss is determined by the radius and material of the cavity, and Q abs is calculated by where Re(β) and |Im(β)| are the real part and the absolute value of the imaginary part of the propagation constant β. For a graphene sheet with a chemical potential of 0.5 eV, when the frequency of the plasmon increases from 80 THz to 83.5 THz, Re(β) and |Im(β)| increase together. However, the rising rate of Re(β) is even much faster than 2|Im(β)|. Consequently, the Q factor increases as the increasing of frequency, and finally results in a similar frequency trend as a function of given g, which are shown in Figure 3a,b. Next, we investigate the frequency splitting ∆F of the corresponding anti-bonding and bonding modes as a function of the azimuthal mode number m of the plasmonic WG mode TM m,1 of the graphene dimers. Figure 5 shows the frequency splitting ∆F of bonding and anti-bonding of WG modes as a function of m, where the distance g is 5 nm. When the azimuthal number increases from 2 to 8, the frequency ∆F of the bonding mode decreases from 0.567 to 0.002 THz, while the counterpart of the anti-bonding mode decreases from 0.275 to 0.0002 THz. When the azimuthal number increases, the wavelength reduces. However, the EM field confinement capacity of the graphene dimer becomes stronger as the resonant wavelength decreases. So the coupling between the plasmonic modes of the nanodisks becomes weaker, which further results in the small frequency splitting for both the bonding and anti-bonding mode, as shown in Figure 5.  Then, we study the frequency and the Q factor of the plasmonic mode TM2,1 of the graphene plasmonic dimer as a function of the chemical potential of the nanodisk at the right side. Figure 6a,b plots the resonant frequency and the corresponding Q factor as a function of μc3 for the four mode categories. In the range from 0.25 to 0.5 eV, one can see that the frequencies of the two modes which are anti-symmetric with respect to the X direction, i.e., [|−1>x,|+1>y] and [|−1>x,|−1>y] modes, are similar. There is a notch in the Q factor of the [|−1>x,|−1>y] mode at μc3 = 0.33 eV, while there is a peak in the counterpart of the [|−1>x,|+1>y] mode at the same point. This is a typical avoid-crossing point [11,14]. However, unlike the photonic dimer composed of the coupled dielectric micro-disks, the Q factors of the [|−1>x,|+1>y] and [|−1>x,|−1>y] modes demonstrate a general climbing trend as the μc3 increases. This is due to the fact that the mismatching between the plasmons inside and outside the nanodisks becomes stronger when μc3 becomes higher, which also prevents the plasmons from escaping from the nanodisk and eventually results in a general increasing of the Q factor. In order to analyze the avoid-crossing property, we study the EM mode evolution as a function of μc3 and the corresponding diagrams are shown in Figure 7. For the [|−1>x,|−1>y] mode, when μc3 reduces from 0.5 to 0.425 eV, the resonant frequency remains stable at a platform of 82.5 THz. Also, the EM field of the left nanodisk does not change significantly and remains in the TM2,1 mode. However, the EM field of the right nanodisk is first weakened, then is enhanced again and switched from the TM2,1 mode to the TM3,1 mode. When μc3 reduces from 0.425 to 0.34 eV, the EM field in the right disk remains stable in the TM3,1 mode. However, the counterpart in the left disk is first weakened then is enhanced again. Also, it switches from the TM2,1 mode to the TM1,1 mode, and the resonant frequency reduces linearly from 83 to 72 THz with μc3. When μc3 further decreases from 0.34 to 0.325 eV, another resonant frequency platform appears at 72 THz. The EM mode of the left disk remains in the TM1,1 mode, and the field in the right disk is first weakened  Then, we study the frequency and the Q factor of the plasmonic mode TM2,1 of the graphene plasmonic dimer as a function of the chemical potential of the nanodisk at the right side. Figure 6a,b plots the resonant frequency and the corresponding Q factor as a function of μc3 for the four mode categories. In the range from 0.25 to 0.5 eV, one can see that the frequencies of the two modes which are anti-symmetric with respect to the X direction, i.e., [|−1>x,|+1>y] and [|−1>x,|−1>y] modes, are similar. There is a notch in the Q factor of the [|−1>x,|−1>y] mode at μc3 = 0.33 eV, while there is a peak in the counterpart of the [|−1>x,|+1>y] mode at the same point. This is a typical avoid-crossing point [11,14]. However, unlike the photonic dimer composed of the coupled dielectric micro-disks, the Q factors of the [|−1>x,|+1>y] and [|−1>x,|−1>y] modes demonstrate a general climbing trend as the μc3 increases. This is due to the fact that the mismatching between the plasmons inside and outside the nanodisks becomes stronger when μc3 becomes higher, which also prevents the plasmons from escaping from the nanodisk and eventually results in a general increasing of the Q factor. In order to analyze the avoid-crossing property, we study the EM mode evolution as a function of μc3 and the corresponding diagrams are shown in Figure 7. For the [|−1>x,|−1>y] mode, when μc3 reduces from 0.5 to 0.425 eV, the resonant frequency remains stable at a platform of 82.5 THz. Also, the EM field of the left nanodisk does not change significantly and remains in the TM2,1 mode. However, the EM field of the right nanodisk is first weakened, then is enhanced again and switched from the TM2,1 mode to the TM3,1 mode. When μc3 reduces from 0.425 to 0.34 eV, the EM field in the right disk remains stable in the TM3,1 mode. However, the counterpart in the left disk is first weakened then is enhanced again. Also, it switches from the TM2,1 mode to the TM1,1 mode, and the resonant frequency reduces linearly from 83 to 72 THz with μc3. When μc3 further decreases from 0.34 to 0.325 eV, another resonant frequency platform appears at 72 THz. The EM mode of the left disk remains in the TM1,1 mode, and the field in the right disk is first weakened Then, we study the frequency and the Q factor of the plasmonic mode TM 2,1 of the graphene plasmonic dimer as a function of the chemical potential of the nanodisk at the right side. Figure 6a,b plots the resonant frequency and the corresponding Q factor as a function of µ c3 for the four mode categories. In the range from 0.25 to 0.5 eV, one can see that the frequencies of the two modes which are anti-symmetric with respect to the X direction, i.e., [|−1> x ,|+1> y ] and [|−1> x ,|−1> y ] modes, are similar. There is a notch in the Q factor of the [|−1> x ,|−1> y ] mode at µ c3 = 0.33 eV, while there is a peak in the counterpart of the [|−1> x ,|+1> y ] mode at the same point. This is a typical avoid-crossing point [11,14]. However, unlike the photonic dimer composed of the coupled dielectric micro-disks, the Q factors of the [|−1> x ,|+1> y ] and [|−1> x ,|−1> y ] modes demonstrate a general climbing trend as the µ c3 increases. This is due to the fact that the mismatching between the plasmons inside and outside the nanodisks becomes stronger when µ c3 becomes higher, which also prevents the plasmons from escaping from the nanodisk and eventually results in a general increasing of the Q factor. In order to analyze the avoid-crossing property, we study the EM mode evolution as a function of µ c3 and the corresponding diagrams are shown in Figure 7. For the [|−1> x ,|−1> y ] mode, when µ c3 reduces from 0.5 to 0.425 eV, the resonant frequency remains stable at a platform of 82.5 THz. Also, the EM field of the left nanodisk does not change significantly and remains in the TM 2,1 mode. However, the EM field of the right nanodisk is first weakened, then is enhanced again and switched from the TM 2,1 mode to the TM 3,1 mode. When µ c3 reduces from 0.425 to 0.34 eV, the EM field in the right disk remains stable in the TM 3,1 mode. However, the counterpart in the left disk is first weakened then is enhanced again. Also, it switches from the TM 2,1 mode to the TM 1,1 mode, and the resonant frequency reduces linearly from 83 to 72 THz with µ c3 . When µ c3 further decreases from 0.34 to 0.325 eV, another resonant frequency platform appears at 72 THz. The EM mode of the left disk remains in the TM 1,1 mode, and the field in the right disk is first weakened then enhanced again; the corresponding mode switches from TM 3,1 to TM 1,2 . When µ c3 further reduces to less than 0.325 eV, the plasmonic field in the right disk remains in the TM 2,1 mode, but the field in the left disk disappears. Meanwhile, the resonant frequency decreases linearly with the reduction of µ c3 . As for the [|−1> x ,|+1> y ] mode, when µ c3 reduces from 0.5 to 0.4 eV, the plasmonic mode in the right disk remains stable in the TM 2,1 mode, while the mode in the left disk is first weakened and then enhanced. Also, the plasmonic mode switches from the TM 2,1 mode to the TM 1,1 mode in the left disk. At the same time, the resonant frequency decreases linearly from 83 THz to 71 THz with the reduction of µ c3 . Similar to the behavior of the [|−1> x ,|+1> y ] mode, a resonant frequency platform appears at around 71 THz when µ c3 decreases from 0.4 to 0.33 eV. In this region, the plasmonic mode of the left disk remains stable in the TM 1,1 mode, and the mode of the right disk is weakened first and then enhanced; also, the TM 2,1 mode is switched to the TM 3,1 mode. When µ c3 further reduces to less than 0.325 eV, the mode field in the left disk weakens gradually and then eventually disappears. However, the right disk remains in the TM 3,1 mode. Meanwhile, the resonant frequency decreases with the reduction of µ c3 . It should be noted that the bonding/anti-bonding state of the dimmer may be changed during the mode evolution procedure [35]. Appl. Sci. 2017, 7, 359 6 of 10 then enhanced again; the corresponding mode switches from TM3,1 to TM1,2. When μc3 further reduces to less than 0.325 eV, the plasmonic field in the right disk remains in the TM2,1 mode, but the field in the left disk disappears. Meanwhile, the resonant frequency decreases linearly with the reduction of μc3. As for the [|−1>x,|+1>y] mode, when μc3 reduces from 0.5 to 0.4 eV, the plasmonic mode in the right disk remains stable in the TM2,1 mode, while the mode in the left disk is first weakened and then enhanced. Also, the plasmonic mode switches from the TM2,1 mode to the TM1,1 mode in the left disk. At the same time, the resonant frequency decreases linearly from 83 THz to 71 THz with the reduction of μc3. Similar to the behavior of the [|−1>x,|+1>y] mode, a resonant frequency platform appears at around 71 THz when μc3 decreases from 0.4 to 0.33 eV. In this region, the plasmonic mode of the left disk remains stable in the TM1,1 mode, and the mode of the right disk is weakened first and then enhanced; also, the TM2,1 mode is switched to the TM3,1 mode. When μc3 further reduces to less than 0.325 eV, the mode field in the left disk weakens gradually and then eventually disappears. However, the right disk remains in the TM3,1 mode. Meanwhile, the resonant frequency decreases with the reduction of μ c3 . It should be noted that the bonding/anti-bonding state of the dimmer may be changed during the mode evolution procedure [35].   As for the modes which are symmetric with respect to the X axis, [|+1>x,|+1>y] and [|+1>x,|−1>y], there is also an avoid-crossing point at around 0.4 eV, which is shown in Figure 6. However, it is not as obvious as that of the asymmetric cases. This may be due to the large frequency splitting, which prevents the strong coupling effect between these two modes. Figure 8 plots the resonance frequency of the [|+1>x,|+1>y] mode as a function of μc3 in a range from 0.25 to 0.5 eV with a gap distance from 0 to 15 nm. The insets show the details in a range between 0.37 and 0.42 eV. Interestingly, with the increasing of the gap distance, the "platform" effect on the resonance frequency becomes more significant, accompanied by the red-shift of the peak position of the Q factor. When the chemical potential μc3 decreases and passes the avoid-crossing point, i.e., 0.4 eV, the EM fields of the two nanodisks react with each other and result in the resonant frequency splitting, which is indicated in Figure 8a, accompanied by the appearance of the "platform" effect, peak of the Q factor, and the switching of the EM field distribution. Nevertheless, when the distance gap g reduces, the coupling becomes stronger, and the "platform" effect becomes weaker, which means that it is not necessary to "accumulate" the frequency reduction to ensure the EM field distribution switching. Also, the peak position of the Q factor shows a blue-shift as shown in Figure 8b.

Conclusions
In this paper, we proposed and numerically analyzed the EM properties of a plasmonic dimer composed of coupled graphene nanodisks. The resonant frequency and the corresponding Q factor of the anti-bonding mode increase as the distance gap increases between the nanodisks, while the opposite trend is observed in the case of the bonding mode. The frequency difference of both the anti-bonding modes and the bonding modes of the TMm,1 reduces as the azimuthal number m increases. An avoid-crossing property has been found in the anti-bonding modes when the chemical potential of the right disk μc3 decreased from 0.5 to 0.25 eV. In the future, more profound graphene nanoclusters will be constructed by adding more nanodisks for the investigation of dynamically tunable Fano resonance lineshape [36] and high order Fano resonance [37], where the graphene nanostructures interact with the external incident EM sources. Another graphene monolayer would be introduced underneath the graphene nanostructures separate ed by a dielectric spacer [38][39][40]. The proposed structure might offer a broad range of applications in both scientific and engineering fields such as light-matter interaction, tunable optical filters, and optical sensors. As for the modes which are symmetric with respect to the X axis, [|+1> x ,|+1> y ] and [|+1> x ,|−1> y ], there is also an avoid-crossing point at around 0.4 eV, which is shown in Figure 6. However, it is not as obvious as that of the asymmetric cases. This may be due to the large frequency splitting, which prevents the strong coupling effect between these two modes. Figure 8 plots the resonance frequency of the [|+1> x ,|+1> y ] mode as a function of µ c3 in a range from 0.25 to 0.5 eV with a gap distance from 0 to 15 nm. The insets show the details in a range between 0.37 and 0.42 eV. Interestingly, with the increasing of the gap distance, the "platform" effect on the resonance frequency becomes more significant, accompanied by the red-shift of the peak position of the Q factor. When the chemical potential µ c3 decreases and passes the avoid-crossing point, i.e., 0.4 eV, the EM fields of the two nanodisks react with each other and result in the resonant frequency splitting, which is indicated in Figure 8a, accompanied by the appearance of the "platform" effect, peak of the Q factor, and the switching of the EM field distribution. Nevertheless, when the distance gap g reduces, the coupling becomes stronger, and the "platform" effect becomes weaker, which means that it is not necessary to "accumulate" the frequency reduction to ensure the EM field distribution switching. Also, the peak position of the Q factor shows a blue-shift as shown in Figure 8b.

Conclusions
In this paper, we proposed and numerically analyzed the EM properties of a plasmonic dimer composed of coupled graphene nanodisks. The resonant frequency and the corresponding Q factor of the anti-bonding mode increase as the distance gap increases between the nanodisks, while the opposite trend is observed in the case of the bonding mode. The frequency difference of both the anti-bonding modes and the bonding modes of the TM m,1 reduces as the azimuthal number m increases. An avoid-crossing property has been found in the anti-bonding modes when the chemical potential of the right disk µ c3 decreased from 0.5 to 0.25 eV. In the future, more profound graphene nanoclusters will be constructed by adding more nanodisks for the investigation of dynamically tunable Fano resonance lineshape [36] and high order Fano resonance [37], where the graphene nanostructures interact with the external incident EM sources. Another graphene monolayer would be introduced underneath the graphene nanostructures separate ed by a dielectric spacer [38][39][40]. The proposed structure might offer a broad range of applications in both scientific and engineering fields such as light-matter interaction, tunable optical filters, and optical sensors.