1. Introduction

Interest in the study of near-field heat transfer has emanated from several theoretical predictions [1,2] and experimental observations [3–6] of unusual behavior in the heat exchanges between object at subwavelength separation distances. These efforts have led to numerous theoretical studies on the effect of geometrical arrangement of objects on the heat transfer [7–10], which have improved our understanding of the far-field and near-field aspects of radiative heat transfer, and revealed the collective effects in many-body systems [11–18]. These newly discovered many-body effects are first of all conceptually very interesting and, by now, they have triggered many theoretical developments on the thermal management in system of objects that exchange heat through radiation. Nevertheless, manipulation of the radiative heat transfer is not simply dictated by the geometrical arrangement of the objects, but also emerge from the precise control of material composition [19–22], shape [23–26], and orientation [14,27–30] of its constituent elements. Exploiting these parameters has led to new ideas in thermal engineering such as the thermal rectification [31–33], localization of thermal modes in fractal structures [13], topological phases in thermal interactions [34–36], thermal barriers [37], modulation of near-field heat transfer using Weyl semimetals [28,38], Magnetic field control of the near-field radiative heat transfer [39], thermal routing [40], and thermal switching [41]. On the basis of theoretical studies, the thermal relaxation process, as well as the thermal steady state distribution, can be tailored by changing the geometric parameters of the system [37,42–44]. A novel methodology, grounded in the linearization of the energy balance equations of the system, demonstrates that the temporal dynamics of temperatures can be simplified to a straightforward eigenvalue problem [45,46]. This approach affords explicit analytical solutions for the temporal progression of temperatures, expressed in the eigenmodes of the system’s response matrix. Recently, this method has been applied to actively control the thermalization process in systems featuring magneto-optical particles [47].

For a comprehensive understanding of both transient and steady-state radiative heat transfer, it is crucial to undertake a thorough investigation into the system’s topological phase. The Su-Schrieffer-Heeger (SSH) model [48], a straightforward lattice configuration showcasing a topological phase transition, has found extensive utility in examining topological characteristics across diverse scientific domains. Despite its simplicity, this model proficiently describes transport properties in both wave and diffusive systems. Our recent investigation into the topological phases of radiative heat transfer has unveiled that a majority of the thermal behaviors in the temporal evolution of temperatures can be attributed to these topological features [49]. These newly uncovered traits in dissipative systems are not only conceptually intriguing but also enable a systematic exploration of the topological phases in terms of the system’s parameters.

In this work, we demonstrate the radiative realization of one-dimensional Su-Schrieffer-Heeger (SSH) model in a finite-size system with an odd number of particles exchanging heat via radiation. The method utilizes the linear response of the system to extract the topological features in terms of the eigenmodes of the response matrix. We show that the topological properties, and consequently the thermal features, can be tailored by changing the volume of particles in a chain of spherical nanoparticles. By confirming the existence of topological phases in the considered system, we investigate the localization of thermal modes as well as the thermalization process in topologically trivial and non-trivial phase of the system.

2. Theoretical model

The radiative analogue of the SSH model we propose in this paper is a chain of nanoparticles with alternating volume, as in Fig. 1(a). Particles have separation distance $d$, with volumes in $A$ and $B$ sub-lattices are given by $V_A=(2-\beta ) V_0$ and $V_B=\beta V_0$, respectively. Here, the number of particles is taken to be odd, $V_0$ is constant, and $\beta \in [0.1,1.9]$ acts as a topological tuning parameter. By these assumptions, the structure will be symmetric with respect to the middle particle ($x_M$), which will be refer to as mirror reflection symmetry. The distance between particles is taken to be constant and large enough compared to the maximum radius $d>4\max \{R_A,R_B\}$. Therefore, the heat transfer could be described in the dipolar approximation regime. The polarizability of the individual nanospheres in the quasistatic approximation is given by $\alpha ^*=3V(\varepsilon -1)/(\varepsilon +2)$. Here, $\varepsilon$ is the dielectric function of the nanospheres. For the polarizability to be consistent with the optical theorem, the radiative damping should be included. Therefore, we use $\alpha =\alpha ^*[1-G_{0}\alpha ^*]^{-1}$ with $G_{0}=ik^3/6\pi$ to address the radiative correction. Here, $k=\omega /c$, and $c$ represents the speed of light. Throughout this paper we will consider Silicon-Carbide (SiC) as typical material, with volumetric heat capacity $c=2.4075\times 10^6~$ JK$^{-1}$m$^{-3}$ and permittivity $\varepsilon =\omega _\infty (\omega _L^2-\omega ^2-i\gamma \omega )/(\omega _T^2-\omega ^2-i\gamma \omega )$, where $\varepsilon _\infty =6.7$, $\omega _T=1.495\times 10^{14}$ rad/s, $\omega _L=1.827\times 10^{14}$ rad/s and $\gamma =0.9\times 10^{12}$ rad/s.

 

Fig. 1. (a) A schematic illustration of a bipartite chain of spherical NPs along the $x$ axis with separation distance $d$, volumes $V_A=(2-\beta ) V_0$ and $V_B=\beta V_0$ immersed in a thermal bath at $T_b=300$ K. (b) Calculated eigenvalue spectrum of chain with $N=61$ NPs, $d=250$ nm, and $V_0=\frac {4}{3}\pi (42)^3$ nm$^3$, for typical values of $\beta \in \{0.5,1,1.5\}$. The eigenvalue spectrum represents a trivial phase ($\beta =0.5$ and $\beta =1$), while it represents topologically nontrivial phase ($\beta =1.5$) with topological edge state inside the topological gap. (c) Eigenvalue spectrum of the radiative SSH model calculated for a chain with $N = 61$ sites, where $d=250$ nm and $V_0=\frac {4}{3}\pi (42)^3$ nm$^3$. Eigenvalues are color-coded based on the Inverse Participation Ratio (IPR) of their corresponding eigenstates, defined as $\texttt {IPR}\mu =\sum _i |\psi \mu (x_i)|^4/\big [\sum _i |\psi _\mu (x_i)|^2\big ]^2$. The bulk energy gap closes at $\beta =1$ and opens again with the appearance of the mid-gap edge states. If $\beta <1$ we are in the trivial phase and the system does not show any topological behaviour at all. On the contrary, if $\beta >1$ we can see the existence of edge states which have topological nature.

Download Full Size | PDF

Now suppose that particles are immersed in a thermal bath at constant temperature $T_b=300$ K. Following previous studies, the differential equation governing the evolution of temperatures is [50,43]

In this equation, $\Theta (\omega,T)$ is the mean energy of a harmonic Planck oscillator and $\tau _{ij}=2{\tt Im}{\mathrm {Tr}}[{\hat {\mathbb {A}}}_{ij}{\tt Im}({\hat {\mathbb \chi }}_{j}){\hat {\mathbb {C}}}_{ij}^{\dagger }]$ represents the radiative energy transmission coefficient between $i$-th and $j$-th NPs [8]. Moreover, 3N $\times$ 3N block matrices are defined as $\hat {\mathbb {A}}=[\hat {\mathbb {1}}-\hat {\alpha }\hat {\mathbb {W}}]^{-1}$, $\hat {\mathbb \chi }=[\hat {\alpha }+\hat {\alpha }G_0^*\hat {\alpha }^{\dagger }]$, and $\hat {\mathbb {C}}=[\hat {\mathbb {W}}+G_0\hat {\mathbb {1}}]\hat {\mathbb {A}}$. Here, $\hat {\mathbb {1}}$ stands for the identity matrix, $\hat {\mathbb {W}}$ is the dipole-dipole interaction matrix and $\hat \alpha$ is 3N $\times$ 3N block diagonal matrix with polarizabilities along the diagonal, i.e., $\hat \alpha = \text {diag}\{ \alpha _A,\alpha _B,\ldots,\alpha _A\}$.

If the difference between the objects temperature and the ambient temperature is small, we can expand $\Theta (\omega, T_j)$ in Eq. (2) about $T_b$ to get

By substituting Eq. (3) into Eq. (1), we get

Now, by defining temperature state vector $\Delta {\textbf T}=({\Delta T}_{1},{\Delta T}_{2},\ldots,{\Delta T_{N}})^T$, we can write Eq. (5) in a compact notation

Here, $C_\mu (0) = \boldsymbol {\phi }_\mu ^\intercal \cdot \Delta {\textbf {T}}(0)$ is the weight of the initial temperature field over $\mu$th mode. One observes that the temporal evolution of temperatures are given in terms of the eigenvalues and eigenstates of the response matrix.

3. Topological phase transition

The eigenvalue spectrum of the response matrix is an essential feature which is independent of the choice of basis. However, it depends sensitively on $\beta$ in the SSH chain of NPs with periodic volumes. To illustrate the topological aspects of the system, the eigenvalue spectrum of the response matrix is shown in Fig. 1(b) for SSH chain of $N = 61$ NPs with lattice constant $d=250$ nm. The spectrum is shown for typical values of $\beta \in \{0.5,1,1.5\}$. Moreover, the volume of NPs are $V_A=(2-\beta ) V_0$ and $V_B=\beta V_0$ with $V_0=\frac {4}{3}\pi (42)^3$ nm$^3$. An important issue to be noted is that, $\beta$ as a coupling parameter, directly influences the eigenvalue spectrum. When $\beta =0.5$, the system can be regarded as a trivial 1D photonic chain with a band gap in eigenvalue spectrum. However, the band gap becomes closed for $\beta =1$. Specifically, when $\beta =1.5$, we observe topological modes with isolated edge as well as bulk states, indicating that the chain is in topologically non-trivial phase.

From the full eigenvalues spectrum in Fig. 1(c), one observes that the system undergoes a topological phase transition from the topologically trivial phase $\beta <1$, where there is no edge state, to the topologically non-trivial phase $\beta >1$ where the edge states appear. It is clear that the spectrum form two bulk bands separated by a gap in both phases and remains gapless only at the critical point $\beta _c=1$. Furthermore, the eigenvalue of the localized edge states appear in the middle of the band gap in the topologically non-trivial phase of the system. It is important to note that the calculations of the eigenvalues and eigenvectors of the response matrix take into account the interaction with all neighboring particles. However, our investigations reveal that the contributions from the nearest neighbors have a dominant influence on both the eigenvectors and the eigenvalue spectrum.

4. Localization and the inverse participation ratio

We observe from Eq. (7) that the eigenstates of the response matrix are good indicator to the topological behaviour of the system. Therefore, we will use the inverse participation ratio (IPR) to measure the degree of localization of these eigenstates. Given an arbitrary state $\psi _\mu (x)$, the IPR is defined as [51]:

The color bars on the right side of Fig. 1(c) depict the IPR of eigenstates in a Su-Schrieffer-Heeger (SSH) chain of nanoparticles (NPs). The eigenvalues displayed in red represents eigenstates with high values of IPR, indicating localization, while the blue color indicates eigenstates with low IPR values, signifying delocalization across the chain. We observe the IPR is maximized ($\sim 0.5$) for the edge states in topologically nontrivial phase, while the eigenstates in bulk bands in both the trivial and nontrivial phases have small IPR, indicating that these modes are mainly extended.

To clarify the localization feature of eigentates in topologically nontrivial phase of the system, we plot in panel (a) of Fig. (2) the eigenstates $\psi _\mu (x)$ where $\mu =1, 2, 3, 4, 31,32,60,$ and $61$ for a typical value of $\beta =1.9$. The profiles are shown as a function of normalized position $x/l$ with $l=x_N-x_1$ for SSH chain with $N=61$ NPs. Since the response matrix is even under parity, the eigenstates are either even or odd with respect to $x/l=0.5$. It is observed that the eigenstate $\psi _1(x)$ is symmetric and spreads uniformly over the entire chain. Since $\psi _1(x)$ has the smallest eigenvalue, it is labeled as $\psi _S(x)$, namely the slowest mode. The eigenstates with larger eigenvalues are alternatively symmetric or asymmetric, and $\psi _{61}(x)$ which has the largest eigenvalue is labeled as $\psi _F(x)$. Further, we observe that the chain supports the existence of degenerate edge states, $\psi _{31}(x)$ and $\psi _{32}(x)$, that are localize at the edges. The edge state located in band with larger [smaller] eigenvalue in topologically trivial phase is referred to as $\psi _U(x)$ [$\psi _D(x)$].

 

Fig. 2. Topological properties for a chain of $N=61$ SiC NPs with $d=250$ nm in topologically nontrivial phase of the system, i.e., $\beta =1.9$. The volume of particles are $V_A=0.1V_0$ and $V_B=1.9V_0$, where $V_0=\frac {4}{3}\pi (42)^3$ nm$^3$. (a) Plot of a bulk states $\{\psi _1,\psi _2,\psi _3,\psi _4,\psi _{60},\psi _{61}\}$ and edge states $\{\psi _{31},\psi _{32}\}$ as a function of normalized position $x/l$. The first column represents symmetric modes and the second column shows the asymmetric modes. (b) The inverse participation ratio (IPR) of eigenstate. (c) Plot of the L-type and R-type eigenstates as a function of normalized position $x/l$.

Download Full Size | PDF

In order to establish the extended or localized nature of the eigenstates, we start by analyzing the IPR of the eigenstates in the topologically non-trivial phase. In Fig. (2(b)), we plot the IPR of all eigenstates when $\beta =1.9$ for $N=61$. It is clearly seen that the IPR is significantly higher ($\sim 0.5$) for two modes $\mu =31$ and $\mu =32$, confirming the existence of a pair of localized edge states, i.e., $\psi _D(x)$ and $\psi _U(x)$, respectively. Instead, the rest of the eigenstates have very small IPR with minimum of $1/61$ for $\mu =1$.

Since $\psi _U(x)$ and $\psi _D(x)$ are the degenerate eigenstates of $\hat H$ in topologically nontrivial phase of the system, we can define hybridized modes $\psi _L(x)=[\psi _D(x)-\psi _U(x)]/\sqrt {2}$ and $\psi _R(x)=[\psi _D(x)+\psi _U(x)]/\sqrt {2}$, namely the L-type and the R-type eigenstates, respectively. As shown in Fig. 2(c), these states are largely concentrated on the left-hand and the right-hand side of the chain, respectively.

We are now in position to explore the influence of the phase transition on the profile of the eigenstates. To this end, we show in Fig. 3(a) the profile of $|\psi _e(x)|^2\equiv |\psi _U(x)|^2 =|\psi _D(x)|^2$ in both the topologically trivial and nontrivial phases across the critical point $\beta _c=1$. It is clear that at $\beta =1.9$, the profile of $|\psi _e(x)|^2$ is completely localized on $x=0$ and $x=l$. However, the component of the edge states on $x=0$ and $x=l$ become smaller as $\beta$ is decreased. As we cross the critical point, these state become extended. We observe that the profile inside the chain is highly uniform in topologically trivial phase and $|\psi _e(x)|^2$ is delocalized over the whole chain. A similar trend is observed for the L-type state, see Fig. 3(b). One sees that the $\psi _L(x)$ is completely localized at $x=x_1$ for $\beta =1.9$. However, this mode delocalizes for smaller $\beta$, and localization collapses as we pass the critical point $\beta _c=1$. As one can see, the state becomes extended for $\beta <1$ and the system does not supports localization feature. Instead, the L-type state has nearly uniform distribution over the entire chain in topologically trivial phase of the system.

 

Fig. 3. Panel (a) Probability distribution $|\psi _e(x)|^2$, and panel (b) Profile of the L-type states in SSH chain of $N=61$ NPs for $\beta \in [0.1,1.9]$, with $d=250$ nm, $V_A=0.1V_0$, $V_B=1.9V_0$, and $V_0=\frac {4}{3}\pi (42)^3$ nm$^3$.

Download Full Size | PDF

5. Classification of topological phases

To achieve a comprehensive understanding of the topological phase transition, this section focuses on examining the topology in reciprocal space. Mathematically representing the response matrix in reciprocal space necessitates defining the configuration of the unit cell. Simillar to the SSH model, the lattice is considered to have periodic boundary conditions. This means that after reaching the last unit cell, the chain continues with the same pattern, creating a cyclic structure. Given that the number of particles is taken to be odd, as illustrated in Fig. 4(a), we can assume a unit cell in which each unit comprises three particles. The particles situated on the faces of each unit cell are shared by two adjacent unit cells. Consider a finite chain and explore two cases, namely $\beta <1$ and $\beta >1$ as illustrated in Fig. 4(a). The key distinction between these configurations lies in $H_{AB}$, which represents the rate of thermal energy transfer from the outermost particle (particles A at both ends) to its nearest neighbor (particles B). It is noteworthy that $|H_{AB}|=|\frac {F_{AB}}{cV_B}|$, indicating its inverse proportionality to volume $V_B$.

 

Fig. 4. (a) Schematic illustration depicting the topological phase transition as the parameter $\beta$ is varied. The system undergoes a significant transformation when $\beta$ transitions from values less than 1 to values greater than 1. For $\beta <1$, particle A (depicted at both ends) efficiently transfers its energy to the chain through particle B due to a substantial energy transfer rate, contributing to the interconnectedness of the system. Conversely, as $\beta$ surpasses 1, particle A becomes decoupled from the chain. The decoupling is visually represented by the isolation of particle A, emphasizing the shift from an interconnected state to a configuration where particle A primarily interacts with the thermal bath. (b) Left: Representation of the hopping terms contributing to the construction of the system’s Hamiltonian. Here, $|H_{ij}|$ represents the energy transfer rate from particle $i$ to particle $j$. Furthermore, $H_{11}$ and $H_{22}$ represent the on-site terms in the Hamiltonian.(b) Right: Representation of inter-cell and intra-cell hopping terms based on the chosen unit cell in the SSH-like Hamiltonian, illustrating the system in reciprocal space. Notably, the specific unit cell selection necessitates the modification (halving) of hopping terms originating from particle A, as depicted in the figure.

Download Full Size | PDF

In the scenario where $\beta < 1$, the volume $V_B$ is comparatively small, leading to a pronounced $|H_{AB}|$. This substantial value of $|H_{AB}|$ facilitates efficient energy transfer from the edge particle to the chain through B. The significant energy transfer rate of the edge particle, driven by this pronounced hopping term, enables an effective channel for energy flow into the chain. Conversely, as the system transitions to $\beta >1$, the edge particle undergoes decoupling from the chain due to the decrease in $|H_{AB}|$. Consequently, it predominantly interacts with the thermal bath. Therefore, the system exhibits localized states at both ends of the chain only for $\beta >1$.

In order to calculate the band structure and study topological phase transitions in SSH-like systems, it is common to use an infinite chain or a chain with a sufficiently large length. The rationale behind this choice is to extract properties that become well-defined and are independent of the specific length of the chain. It’s important to clarify that, for the purposes of our analysis, we employ the infinite chain as a theoretical approximation to facilitate the extraction of the band structure and the topological characteristics of the system. The second-quantized Hamiltonian of the system can be represented in terms of energy transfer hoppings. In this framework, the hopping parameters characterize the energy transfer rates between particles or between particles and the thermal bath. As illustrated in Fig. 4(b), the terms $H_{ij}\equiv c^{-1}F_{ij}V_j^{-1}$, represent the hopping coefficients from particle $i$ to particle $j$ in the system. Moreover, given that $F_{ij}=F_{ji}$, it can be demonstrated that $H_{ij}=H_{ji}$ holds for even values of $|i-j|$. However, $H_{ij}\neq H_{ji}$ if $|i-j|$ is odd. The hoppings with $|i-j|=0$ correspond to the cooling power contributions from the respective objects, influencing the on-site terms of the Hamiltonian, akin to the Rice-Mele model [52]. On the other hand, the hoppings with $|i-j|=1,2,\ldots$ represent the hopping between first, second, $\cdots$ neighbors, respectively. It should be emphasized that the influence of $H_{ij}$ with $|i-j|>1$ on the band structure of the Hamiltonian is deemed negligible. These terms are significantly smaller than the nearest-neighbor hoppings. However, neglecting these terms does not affect the consideration of many-body interactions between particles, as this effect is included in the hopping terms through $F_{ij}$ in Eq. (3).

Considering our unit cell configuration, the particles located on the faces of each unit cell are shared by two adjacent unit cells. Therefore, if we intend to describe the inter-cell and intra-cell interactions of our system, the half hopping should be assigned to the transmission probabilities involving these particles. The hopping parameters for this specific choice of unit-cell are depicted on the right-hand side of Fig. 4(b), where all the hoppings originating from particle A are halved. The second-quantized Hamiltonian up to the 4th neighbor is then given by:

Here, $a_n$, $b_n$ and $c_n$ are the annihilation operators on sites A, B and A in a unitcell, respectively, in the n-th unit cell. The matrix representation of $\mathcal {\hat H}$ is given by:

 

Fig. 5. (a) Band structure of the model. The plot illustrates the eigenvalue spectrum as a function of $k$ in the first Brillouin zone. The $+$ and $-$ symbols represent the upper and lower energy bands, respectively. A band gap is clearly visible for $\beta \neq 1$, signifying the presence of distinct energy bands. Notably, the band gap is closed when $\beta =1$, where the upper and lower bands touch, marking a critical point in the topological phase transition. (b) Zak phase analysis of the model. The plot depicts the Zak phase as a function of the parameter $\beta$, showcasing the topological phase transition in the system. For $\beta <1$, both upper and lower bands exhibit a Zak phase of zero. As $\beta$ crosses the critical value of 1, a dramatic change occurs, leading to Zak phases of $+\pi$ for the upper band and $-\pi$ for the lower band.

Download Full Size | PDF

The eigenstates of the upper and lower bands satisfy $\mathcal {\hat H} (k)|\psi ^\pm (k)\rangle =\lambda _k^\pm |\psi ^\pm (k)\rangle$. Correspondingly, we denote the left eigenstates as $\langle \phi ^\pm (k)|$ with $\langle \phi ^\pm (k)|\mathcal {\hat H} (k)=\langle \phi ^\pm (k)|\lambda _k^\pm$. It is noteworthy that these eigenvectors play a crucial role in the investigation of the topological phase transition, particularly through the analysis of the Zak phase of the system. The Zak phase, an extension of the Berry phase [53], is a fundamental concept applicable to Bloch wave functions within solid materials. Specifically, it quantifies the phase change of the periodic part of the Bloch wave function as the wave vector $k$ traverses the Brillouin zone. In essence, the Zak phase provides a geometric characterization of the global properties of the Bloch states in reciprocal space, offering valuable insights into the topological nature of the system. The geometric Zak phase can be determined for both the upper and lower bands from

As depicted in Fig. 5(b), we explore the Zak phase behavior for both bands concerning the parameter $\beta$. In the regime where $\beta < 1$, the Zak phases for both bands are consistently identified as zero. This observation suggests that, in the $\beta < 1$ regime, the system resides in a topologically trivial phase. However, as we traverse the parameter space and surpass $\beta > 1$, a notable topological phase transition unfolds. At this point, the Zak phases for the upper and lower bands become $+\pi$ and $-\pi$, respectively. This substantial change in Zak phases serves as a distinctive indicator of the emergence of a topologically non-trivial phase within the system. The transition is closely linked to the development of edge states localized at the boundaries of the system.

6. Robustness of edge states

The robustness of topological edge states refers to their ability to remain immune to disorder or imperfections in the system and maintain their protected nature. As illustrated in Fig. 1(c), the presence of a gap in the eigenvalue spectrum distinguishes the edge states from the bulk states. This gap serves as a protective barrier, preventing hybridization between the edge states and the bulk states. Consequently, the edge states exhibit resilience against local perturbations that may lead to the closure of the gap. Simultaneously, there is an expectation that the edge states remain symmetric/asymmetric, allowing us to construct the L-type and R-type eigenstates that exhibit localization properties within a topologically robust structure. Hence, the robustness of a topological edge state is affirmed when two essential criteria are met. Firstly, the gap between the edge states and the bulk states remains unaltered in the presence of external perturbations. Secondly, the preservation of the edge state’s degeneracy under the influence of perturbations is crucial, indicating that the degeneracy remains intact and unaffected by external factors. It is noteworthy that decrease in the IPR of the edge states due to perturbation may suggest a reduction in their robustness. However, it is important to clarify that changes in the IPR of edge states do not necessarily imply a lack of robustness. Rather, any change in the IPR of a state indicates a shift in the degree of its localization. These criteria serve as indicators of the robustness of the topological edge states employed in this section.

6.1 Response matrix perturbation

In order to provide a general evidence for the robustness of topological edge states in the proposed system, we investigate the influence of diagonal, off-diagonal and symmetric disorder on the eigenvalue spectrum of the response matrix. The disorder is induced by substituting a target element in the response matrix, denoted as $H_{ij}$, with a random variable sampled from a Gaussian distribution $f(h)=\exp [-(h-H_{ij})^2/2\sigma ^2]/\sqrt {2\pi }\sigma$. Therefore, denoted concisely as $H_{ij}\to h\sim \mathcal {N}(H_{ij},\sigma ^2)$, the parameter $\sigma$ serves as a quantification of the magnitude of local disorder, encompassing factors such as displacement, removal, and changes in material or particle volume. The results presented in this section are obtained for $\beta =1.9$ and are averaged over 500 realizations to ensure statistical significance and reliability.

The diagonal perturbation for each sample is implemented by applying a Gaussian distribution $\mathcal {N}(H_{ii},\sigma ^2)$ to all the diagonal elements of the response matrix. As shown in Fig. 6(a), introducing diagonal disorder disrupts the sublattice symmetry. This leads to a loss of degeneracy in the edge states and the closure of the relaxation gap. Consequently, this leads to the loss of topological protection for the edge modes within the system.

 

Fig. 6. Effect of (a) diagonal, (b) off-diagonal, and (c) mirror reflection symmetric disorder on the eigenvalue spectrum of the response matrix in a finite chain comprising $N = 61$ particles. The eigenvalue spectra are presented as a function of the disorder strength parameter $\sigma$. Gaussian perturbations are applied to the elements of the matrix to introduce disorder, and the presented results are averaged over 500 realizations for each value of the disorder strength parameter $\sigma$.

Download Full Size | PDF

Figure 6(b) illustrates the impact of off-diagonal disorder on the eigenvalue spectrum. In this case, perturbation for each sample is applied by subjecting all the off-diagonal elements of the response matrix to a Gaussian distribution $\mathcal {N}(H_{ij},\sigma ^2)$. Although the degeneracy of the edge states is preserved in this case, the gap narrows with an increase in the disorder strength parameter. Based on the results, the localized edge states exhibit resilience to this specific perturbation for small perturbation strengths, specifically when $\sigma \lesssim 50$. In both scenarios, the disorder parameter $\sigma$ exerts a notable influence on the eigenvalues of the bulk bands. As the magnitude of the disorder approaches the band gap, the gap completely closes.

In Fig. 6(c), we present the results of mirror reflection symmetric disorder applied to the elements of the response matrix. In this scenario, the perturbation is applied by subjecting the elements of the response matrix, denoted as $H_{ij}$ and $H_{ji}$ with $i,j\leq N/2$, to a Gaussian distribution. Subsequently, the same values are used for the elements $H_{N-i,N-j}$ and $H_{N-j,N-i}$, respectively, to ensure the preservation of mirror-reflection symmetry. Due to the preservation of this symmetry, which is inherent to the underlying topological properties, we observe the persistent presence of both the gap and degeneracy in the eigenvalue spectrum even under substantial perturbation magnitudes.

6.2 Disorder induced by perturbations in positions or volumes

In the previous section, we applied uncorrelated perturbations to the elements of the response matrix. However, in practical scenarios, disorders like positional displacements or changes in volumes can affect all elements of the response matrix. This can result in correlated, non-Gaussian perturbations in these elements. To address this phenomenon, the current section focuses on specific categories of disorder, such as perturbations in particle positions or volumes, with the aim of gaining a comprehensive understanding of their impact on the robustness of the edge states. The results presented in this section are derived with $\beta =1.5$ and represent the average of 200 realizations in each case, ensuring both statistical significance and reliability.

6.2.1 Perturbations in positions

The perturbations arising from the disorder in the spatial arrangement of particles within the system are induced by replacing the positions of particles, denoted as $x_i$, with randomly sampled variables from a Gaussian distribution $f(x)=\exp [-(x-x_i)^2/2\sigma _x^2]/\sqrt {2\pi }\sigma _x$. Consequently, the parameter $\sigma _x$ quantifies the magnitude of position disorder, providing a meaningful measure of its strength. This parameter signifies the extent to which the positions of particles deviate from their nominal values due to the introduced perturbation.

The impact of disorder in the particle positions on the eigenvalues of the eigenstates is depicted in Fig. 7(a). As this disruption breaks the underlying topological invariant, namely Mirror-reflection symmetry, it therefore disrupts the symmetric/asymmetric nature of the eigenstates. The edge states are no exception to this rule, and as a consequence of this perturbation, the degeneracy in their eigenvalues is lifted. Therefore, the edge states will not exhibit robustness against this type of perturbation. As depicted in Fig. 7(b), it is observed that the IPR of the two edge states diverges as the perturbation strength is increased. Therefore, the perturbation instigates a decrease in the confinement of the edge state, signifying the transition of the edge states towards delocalization within the bulk of the material. Additionally, the perturbation amplifies the IPR of the bulk states, signifying the emergence of localized states within the bulk material.

 

Fig. 7. Effect of Gaussian perturbations on the eigenvalue spectrum and inverse participation ratio (IPR) in a 61-particle chain with $\beta =1.5$. Gaussian disorder is applied uniformly to all particles in each case. Results are averaged over 200 realizations for each $\sigma$. Panels (a), (b) depict perturbations with a Gaussian distribution on particle positions. Panels (c), (d) showcase mirror-reflection symmetric perturbations on particle positions. Panels (e), (f) illustrate perturbations with a Gaussian distribution on particle volumes. Panels (g), (h) present mirror-reflection symmetric perturbations on particle volumes.

Download Full Size | PDF

Mirror-reflection symmetric disorder is implemented by subjecting the positions of particles in the left half of the chain to a Gaussian distribution $\mathcal {N}(x_i,\sigma _x^2)$, and applying a mirror replacement to the other half. As shown in Fig. 7(c),the edge state degeneracy remains unaffected by this perturbation for $\sigma _x\lesssim 10$, demonstrating the robustness of the topological edge state. However, as expected, for large perturbations, the degeneracy unfolds, and the gap closes. It is important to note that the perturbation also influences other modes of the system, as in the previous case. As demonstrated in Fig. 7(d), the perturbation strength’s influence on the IPR of the eigenstates of the response matrix is shown. From this figure, it is evident that the perturbation does not have a significant effect on the localization of the edge states for $\sigma _x\lesssim 10$. On the contrary, the disruption of translational symmetry results in the emergence of additional localized sites within the chain. An increase in the IPR of bulk modes serves as an indicative sign of the formation of these localized sites inside the chain.

6.2.2 Perturbations in volumes

The robustness of edge states against perturbations in particle volumes is investigated by introducing disorder through the substitution of the volumes of particles, denoted as $V_i$, with randomly sampled variables from a Gaussian distribution $f(V)=\exp [-(V-V_i)^2/2\sigma _v^2]/\sqrt {2\pi }\sigma _v$. Consequently, the parameter $\sigma _v$ quantifies the extent to which the volumes of particles deviate from their nominal values. Furthermore, it is crucial to note that, akin to the position perturbations, this perturbation has the potential to lead to correlated non-Gaussian disturbances within the response matrix elements. In Figure 7(e), it is evident that the system lacks topological robustness against this specific perturbation. While this perturbation preserves the inverse participation ratio (IPR) of the edge states for small perturbation strengths (i.e., $\sigma _v\lesssim 0.2 V_0$ as shown in Fig. 7(f)), it results in the unfolding of the degeneracy of the edge states.

The influence of mirror-reflection symmetric perturbations in volumes on the eigenvalue spectrum and inverse participation ratios (IPRs) of the modes is depicted in Fig. 7(g) and Fig. 7(h), respectively. We observe that the edge states demonstrate robustness against this type of perturbation. Concurrently, we observe an increase in the localization of the remaining states and an eventual closure of the gap with an increase in perturbation strength, consistent with expectations.

7. Thermalization process

We now investigate how the choice of the coupling parameter $\beta$ influences the system’s thermalization process. To this end, the transient regime of the solutions for $\Delta T(x_1,t)$ is shown in Fig. 8(a). The initial condition for temperatures are defined such that, only particle $1$ is heated up to $350$ K, while all other particles being initially at room temperature $300$ K, i.e., $\Delta T(x,0)=50\delta (x-x_1)$. We use the value of $\beta =0.5$ and $\beta =1.9$, which will give topologically trivial and non-trivial phases, respectively. Clearly, the thermalization process depends on the topological phase of the system. A close comparison of the results indicates that the temperature $\Delta T_1(t)$ decreases much faster in the nontrivial phase (solid-blue), than that for the trivial phase (dashed-green). In the topologically non-trivial phase, particle $1$ can be maximally coupled with the thermal bath with isolated relaxation dynamics. In contrary, the thermalization process in the trivial phase of the system could be very slow due to the contribution of the extended states. The ultrafast thermalization dynamics observed in Fig. 8(a) can be understood by considering the temporal evolution of temperatures for the initial condition $\Delta T(x,0)=50\psi _L(x)$ in topologically nontrivial phase of the system. Since $\psi _L(x)\sim \delta (x-x_1)$ in this phase, we observe a similar and ultrafast thermalization dynamics for excitation of the L-type eigen state.

 

Fig. 8. (a) Thermalization process in topologically trivial ($\beta =0.5$) and nontrivial ($\beta =1.9$) phase of the system. The blue (dashed-green) curve represents the transient evolution of $\Delta T(x_1,t)$ in topologically nontrivial (trivial) phase of chain with $N=61$ nanoparticles subjected to an initial temperature condition $\Delta T(x,0)=50\delta (x-x_1)$. The red curve shows the transient evolution of $\Delta T(x_1,t)$ under initial condition $\Delta T(x,0)=50\psi _L(x)$ in topologically nontrivial phase of the system. (b) The weight distribution of the initial temperature vector over modes for each case.

Download Full Size | PDF

In order to understand the physics behind the thermalization process, the weight of the initial temperature field are shown in Fig. 8(b), as a function of $\mu$. We find that heating up particle $1$, can sufficiently excite the edge states in the topologically nontrivial phase. It is clear that the weight distribution of the initial temperature state is localized at $\lambda _e$, similar to the $\psi _L(x)$. On the other hand, heating up particle $1$ in non-trivial phase, mainly excites the localized L-type state. Consequently, similar to the L-type initial temperature state, the temperature $\Delta T_1$ decays as $\Delta T_1(t)\sim 50\exp (-\lambda _e t)$. In contrary, heating up particle $1$ in the topologically trivial phase, excites most of the modes in the bulk bands. Therefore, the contribution of the slowest mode $\mu =1$ results in slow thermalization process in this case.

To compare the thermalization process between topologically trivial and non-trivial phases, in panel (a) of Fig. 9, we present a detailed comparison of the relaxation process for $\beta =1.9$ (non-trivial phase) and $\beta =0.1$ (trivial phase), where only the particle at the edge is initially heated to $\Delta T_1(0)=50$ K. In the non-trivial phase, particle 1 is relatively small, allowing for its rapid relaxation with minimal interaction with the chain. However, in the trivial phase, the interaction of particle 1 with the chain through particle 2 is significant due to the smaller coupling of particle 1 with the thermal bath at this phase. This competition between interaction with the bath and interaction with the chain results in a noticeable variation in the temperatures of the other particles during the thermalization process, as illustrated in the figure.

 

Fig. 9. Comparison of thermalization processes in topologically trivial and non-trivial phases in a chain with 61 nanoparticles. Panel (a) illustrates the relaxation process for $\beta =1.9$ (non-trivial) and $\beta =0.1$ (trivial) phases, initiated by heating only the edge particle to $\Delta T_1(0)=50$ K. In the non-trivial phase, the small size of particle 1 leads to rapid relaxation with minimal chain interaction. However, in the trivial phase, significant interaction occurs between particle 1 and the chain, resulting in notable temperature variations in other particles. Panel (b) explores a different scenario by initiating the heating process with an increase in temperature of the second chain particle, $\Delta T_2(0)=50$ K. In the non-trivial phase, the large size of particle 2 leads to prolonged and substantial chain interaction, affecting the temperature of the edge particle. In contrast, the rapid cooling of this particle, attributed to its smaller volume in the topologically trivial phase, inhibits its interaction with the chain. Panel (c) Thermalization process in both topologically trivial and non-trivial phases, with the edge particle’s temperature fixed at $\Delta T_1 = 50$ K. The energy flux through the chain results in an increase in temperature, with smaller equilibrium values of temperature within the chain observed in the topologically non-trivial phase. Panel (d) Evolution of temperatures toward thermal equilibrium with the second particle’s temperature fixed at $\Delta T_2 = 50$ K. Large equilibrium temperatures are observed in both phases, indicating a higher heat flux within the chain in both scenarios.

Download Full Size | PDF

Moving to panel (b) of Fig. 9, we replicate the calculations with a different scenario. This time, we initiate the heating process by increasing the temperature of the second particle in the chain, denoted as $\Delta T_2(0)=50$ K. In the topologically non-trivial phase, the large size of this particle leads to reduced radiative cooling, resulting in a more prolonged and substantial interaction with the rest of the particles in the chain. Consequently, we observe that the temperature of particles in the chain, especially the edge particle (i.e., $\Delta T_1$), is also significantly affected during the thermalization process.

For a more comprehensive understanding of the influence of topological phases on the thermal properties, we repeated the previous scenarios with a variation. This time, the temperatures of the targeted particles are held fixed during the thermalization process, intending to extend the interaction time and observe the distinctive effects of topology. In panel (c) of Fig. 9, the thermalization process is compared between the two phases, with the temperature of the particle at the edge fixed at $\Delta T_1=50$ K. In both phases, the energy flux through the chain leads to an increase in temperature. However, it is noteworthy that the maximum equilibrium values of temperature within the chain ($\Delta T_2$) differ between the two phases. These results suggest a more pronounced energy transfer and temperature increase inside the chain in the topologically trivial phase compared to the non-trivial phase. This observation is consistent with the findings in panel (a).

The results for the temperature evolution toward thermal equilibrium, with the temperature of particle 2 fixed at $\Delta T_2=50$ K, are presented in panel (d) of Fig. 9. In this scenario, we observe that the equilibrium temperatures are higher compared to the previous case. In the non-trivial phase, the large volume of particle 2 facilitates significant radiative transfer through the chain, resulting in a large steady-state temperature $\Delta T_1(t\to \infty )\sim 50$ K for the edge particle. Conversely, in the topologically trivial phase, maintaining the temperature $\Delta T_2$ fixed allows us to observe the flow of energy in the system, in contrast to the relaxation process observed in Fig. 9(b). It’s noteworthy that in both scenarios, a substantial heat flux occurs within the chain, as observed in panel (c) for the topologically trivial phase.

8. Conclusion

In summary, our study elucidates the presence of topological modes encompassing both edge and bulk states, along with a discernible topological phase transition in a chain featuring an odd number of particles exchanging energy via radiation. Commencing from the energy balance equation, we harnessed the response matrix as a pivotal tool to delineate the emergent topological phases. Notably, we have established that within the topologically nontrivial phase, the thermalization process is expedited for edge particles, whereas it proceeds at a more measured pace for particles situated within the bulk. Furthermore, the influence of perturbations underscores the resilience of the edge state in the face of disorders characterized by mirror-reflection symmetry. The robustness of our topological analysis approach showcased in this work may serve as a catalyst for future investigations in the realm of radiative heat transfer, offering a promising avenue for the exploration of thermal topology in the pursuit of steadfast heat transfer mechanisms.