Dynamical classification of topological quantum phases
Lin Zhang Long Zhang Sen Niu Xiong-Jun Liu




d
i
𝕀1





These authors contribute equally to this work.
International Center for Quantum Materials, School of Physics, Peking University, Beijing 100871, China.
Collaborative Innovation Center of Quantum Matter, Beijing 100871, China.


These authors contribute equally to this work.
International Center for Quantum Materials, School of Physics, Peking University, Beijing 100871, China.
Collaborative Innovation Center of Quantum Matter, Beijing 100871, China.


International Center for Quantum Materials, School of Physics, Peking University, Beijing 100871, China.
Collaborative Innovation Center of Quantum Matter, Beijing 100871, China.


International Center for Quantum Materials, School of Physics, Peking University, Beijing 100871, China.
Collaborative Innovation Center of Quantum Matter, Beijing 100871, China.
Institute for Quantum Science and Engineering and Department of Physics, Southern University of Science and Technology, Shenzhen 518055, China


Topological phase of matter is now a mainstream of research in condensed matter physics, of which the classification, synthesis, and detection of topological states have brought excitements over the recent decade while remain incomplete with ongoing challenges in both theory and experiment. Here we propose to establish a universal dynamical characterization of the topological quantum phases classified by integers, and further propose the high-precision dynamical schemes to detect such states. The framework of the dynamical classification theory consists of basic theorems. First, we uncover that classifying a d-dimensional (dD) gapped topological phase of generic multibands can reduce to a (d-1)D invariant defined on so-called band inversion surfaces (BISs), rendering a bulk-surface duality which simplifies the topological characterization. Further, we show in quenching across phase boundary the (pseudo)spin dynamics to exhibit unique topological patterns on BISs, which are attributed to the post-quench bulk topology and manifest a dynamical bulk-surface correspondence. For this the topological phase is classified by a dynamical topological invariant measured from dynamical spin-texture field on the BISs. Applications to quenching experiments on feasible models are proposed and studied, demonstrating the new experimental strategies to detect topological phases with high feasibility. This work opens a broad new direction to classify and detect topological phases by quantum dynamics.



Dynamical classification of topological quantum phases
    Xiong-Jun Liu
[Correspondence author: xiongjunliu@pku.edu.cn]
    December 30, 2023
=================================================================



In the last over thirty years, following the discovery of the quantum Hall effect <cit.>, physicists have established a new classification of fundamental states of quantum matter, called topological quantum phases <cit.>, which are beyond the scope of the celebrated Landau-Ginzburg-Wilson framework that characterizes quantum matter by local symmetry-breaking orders <cit.>. A topological quantum phase can be classified by bulk topological invariant defined in the ground state at equilibrium, and hosts protected boundary modes through the bulk-boundary correspondence <cit.>.
This characterization crucially influences the strategies of identifying topological states in experiment, which is challenging in general since a non-local topological invariant may not have direct physical measurements.
Despite the great success achieved with various strategies in discovering new topological matter, such as topological insulators <cit.> and semimetals <cit.> whose boundary modes can be probed by transport measurement or angle resolved photoemission spectroscopy, by far only a small portion of topological states predicted in theory have been uncovered in experiment <cit.>. In some circumstances the measurements are not fully unambiguous for being not direct observations of topological numbers, for which there are topological states, e.g., the topological superconductors <cit.>, even well explored in theory, necessitating further experimental confirmation.

As yet, the topological states have been primarily characterized with classification theories developed for equilibrium systems, ranging from symmetry-protected topological phases <cit.> to intrinsic topological orders <cit.>.
A fundamental question arises that, for a generic topological quantum phase defined for the equilibrium ground state, is there a non-equilibrium classification of such phase? We address this fundamental issue in the present work, and establish a dynamical classification theory for topological quantum phases characterized by integers of multiband and all dimensions. Not only being of the clear fundamental significance, the dynamical classification is also experimentally important in probing new topological quantum physics, particularly relevant for ultracold atom systems where, due to heating, the equilibrium ground states are usually hard to achieve, but the quantum dynamics can be readily engineered as having been considered in the recent studies of topological systems <cit.>. The cold atom experiments have demonstrated the feasibility of simulating novel topological systems with quantum gases, including the 1D Su-Schrieffer-Heeger (SSH) chain <cit.>, 1D chiral topological phase <cit.>, and 2D Chern insulators <cit.>. In particular, characterizing the band topology of an equilibrium Hamiltonian from quantum dynamics was proposed and explored very recently in theory and experiment <cit.>. Nevertheless, these studies are not universal and valid for two-band models in specific dimensions. No dynamical classification theory was proposed until the present work.

The framework of the dynamical classification theory is made up of several basic theorems uncovered in this work. First, for a generic d-dimensional (dD) gapped topological phase characterized by integer invariants, we show a bulk-surface duality that its classification reduces to a (d-1)D invariant defined on the so-called band inversion surfaces (BISs). Further, in quenching from trivial to topological phases, the BISs are captured dynamically that on BISs the time-averaged (pseudo)spin-polarizations vanish in the unitary evolution. Finally, the bulk topology of the post-quench phase is classified by the dynamical topological invariant determined via an emergent dynamical spin-texture field on BISs. The results manifest a dynamical bulk-surface correspondence which can be directly measured with high feasibility in quenching experiments. With the dynamical classification theory we propose various high-precision dynamical schemes to detect topological quantum states, of which the application to measuring Chern insulators has been remarkably achieved in a latest experiment <cit.>.

Generic model.–We start with the generic d-dimensional (dD) gapped topological phases, including insulators and superconductors, which are classified by integer invariants in the Altland-Zirnbauer (AZ) symmetry classes <cit.>. The basic Hamiltonian can be written in the elementary representation matrices of the Clifford algebra <cit.> (see also Supplementary Material <cit.>)

    ℋ(𝐤)=h⃗(𝐤)·γ⃗=∑_i=0^dh_i(𝐤)γ_i,

where the γ⃗ matrices define a (pseudo)spin obeying anti-commutation relation {γ_i,γ_j}=2δ_ij𝕀, with i,j=0,1,…,d, and h⃗(𝐤) mimics a (d+1)D Zeeman field depending on the Bloch momentum k in BZ.
The dimensionality of γ matrices in the elementary representation reads n_d=2^d/2 (or 2^(d+1)/2) if d is even (or odd), which is the minimal requirement to open a topological gap for the dD topological phase <cit.>. In the 1D/2D regimes, for instance, the γ matrices simply reduce to the Pauli ones, and ℋ(𝐤) describes a two-band model for the topological states, e.g., the well-known SSH model <cit.> for 1D BDI class insulator and the Haldane model <cit.> for 2D Chern insulator. Similarly, for 3D/4D phases, the γ matrices take the Dirac forms, and a fully gapped topological phase has to involve at least four bands, such as the 3D AIII class topological insulator (also DIII class supercondoctors) <cit.> and 4D quantum Hall effect <cit.>. While for convenience the dynamical classification theory is formulated with the basic Hamiltonian written in the above elementary representation, the theory applies to any generic multiband model of the dD phase (see proof in Appendix). Therefore, all the main results presented based on Eq. (<ref>) directly apply to the generic dD gapped topological phases characterized by integer invariants.




To facilitate the description, we pick up an arbitrary component, say h_0(𝐤), from the vector field h⃗(𝐤) to characterize the `dispersion' of the n_d decoupled bands. Accordingly, we denote the remaining components of h⃗(𝐤) as the `spin-orbit' (SO) field h_ so( k)=(h_1,...,h_d). Without the SO field, the band crossing occurs for the n_d decoupled bands if h_0( k)=0 results at certain momentum points of the BZ. All such momentum points form closed (d-1)D surface, dubbed band inversion surface (BIS), which is a key concept in this work. The BIS is related to, but conceptually different from the
familiar phenomenon of band inversion in topological insulators, and can be defined for broad classes of phases including topological superconductors. Across the BIS the energy difference between half of the n_d bands and the remaining half switches sign. In general one can have multiple BISs in the first BZ. The gap opens when the nonzero SO field h_ so( k) is switched on in the BISs. As shown below, the dD gapped phase being topologically nontrivial necessitates that the BZ includes at least one BIS which hosts a nonzero (d-1)D invariant.

Bulk-surface duality.–The topology of the Hamiltonian H( k) is classified by dD winding number for odd dimensionality d=2n-1, or the n-th Chern number if d=2n. The topological characterization of the equilibrium phase is formulated by a mapping from the BZ, i.e. a dD torus T^d, to the dD spherical surface S^d through the unit vector field n( k)=h⃗/|h⃗|. The topological number counts the times that the mapping covers S^d. Interestingly, as detailed in the Appendix, we show that the topological characterization by dD winding number (or Chern number) reduces to a (d-1)D Chern number (or winding number) w_d-1 defined on the BISs

    ν_d(h_0, h_ so)=w_d-1( h_ so),

where w_d-1 characterizes the integer times that the unit SO field ĥ_ so( k)=h_ so/|h_ so| covers over the (d-1)D spherical surface S^d-1 when k runs over all the BISs of the BZ. Without giving details, we show the above invariant w_d-1 intuitively by clarifying the key ingredients of the theorem, while the rigorous proof can be found in the Appendix. As for topological phase, the mapping T^d↦ S^d covers nonzero integer times of the dD spherical surface, any component of the vector field h⃗( k), including h_0( k), spans both positive and negative values in the BZ. Thus the BIS with h_0( k)=0 must be obtained. It is convenient to denote by V_ BIS the vector volume of the BZ with h_0( k)<0, and V̅_ BIS the volume with h_0( k)>0. Note that w_d-1 in the right hand side of Eq. (<ref>) is given by a (d-1)D surface integral enclosing the open vector volume V_ BIS (Appendix), which is gauge-independent and counts the total topological charges located at the singularities with h_ so( k)=0 in this volume (see Fig. 1a). This result can be intuitively understood in the following way, for which we rewrite that h_0( k)=m_0+h̃_0( k), with h̃_0 being momentum-dependent and m_0 mimicking a constant `magnetization' which tunes the topology of the phase. It is ready to know that when m_0>max[|h̃_0( k)|] for all k, the phase is trivial since no BIS can be obtained in the BZ. The BIS starts to emerge in the BZ by reducing the magnetization to m_0≲max[|h̃_0( k)|, while the gap of the system keeps open if no topological charges (at h_ so( k)=0) enter the vector volume V_ BIS. Further reducing m_0 eventually enables that the topological charges pass through the BIS and enter the volume V_ BIS. Note that when a topological charge crosses the BIS, the bulk gap closes one time since h_0= h_ so=0 occurs at the crossing momentum on the BIS. This leads to a topological phase transition with the topological invariant varying by one. Then the number of topological charges with h_ so( k)=0 enclosed in the volume V_ BIS reflects the times of gap closing by tuning m_0 from trivial to the final topological regime, and the total topological charge gives the bulk topological number of the phase, as characterized by the theorem in Eq. (<ref>).

The above theorem shows a bulk-surface duality which maps the classification of bulk topology to the characterization on BISs. The bulk-surface duality provides a simplified characterization of the topological phases with lower-dimensional invariants, while there are two important issues having to be addressed. First, the BISs are not unique but depends on the choice of h_0 axis. Then how to precisely identify BISs in a real experiment? Further, how to measure the topology defined on the BISs? As we show below that, the uncovered bulk-surface duality is further mapped to a novel dynamical form in the quench dynamics, which directly measures all the topological features on BISs with high feasibility.


Non-equilibrium classification: dynamical bulk-surface correspondence.–We proceed to consider the non-equilibrium characterization of topological phases based on the bulk-surface duality shown above. The quantum dynamics can be induced by quenching the phase from a deep initially trivial phase at time t=0, which is given by setting m_0|_t<0≫max[|h⃗( k)|], to the final topological state which is obtained in proper regime |m_0||_t>0<max[|h_0( k)|] (Fig. 1b). The unitary evolution of the (pseudo)spin γ⃗ is then governed by the post-quench Hamiltonian H(m_0, k). The quenching dynamics can be quantified by introducing the dynamical averaging of the spin-polarization, which is obtained by

    ⟨γ_i( k)⟩=lim_T→∞1/T∫_0^T t Tr[ρ_0 e^ℋtγ_ie^-ℋt],

where ρ_0 is density of matrix of the initial state. The peculiar effects are obtained on the BISs. Note that the system is initialized in the fully spin-polarized phase with the (pseudo)spin of all Bloch states pointing along negative h_0 axis (for positive m_0). The spin processes with respect to h⃗( k) after quenching. On the BISs, the vector h⃗( k)= h_ so( k) is perpendicular to initial spin polarization, which leads to the procession of the spin vector γ⃗ within the plane perpendicular to h_ so while incorporating h_0( k) axis (Fig. 1c). As a result, the dynamical averaging ⟨γ_i⟩ vanishes right on the BISs for all i-components. On the other hand, at the momentum k away from the BISs the h_0 component is nonzero, so h⃗ is not perpendicular to the initial spin polarization. Then the spin procession leads to a nonzero ⟨γ⃗⟩ along the h⃗( k)-direction (Fig. 1c). With these results we conclude a dynamical characterization of the BISs that

    ⟨γ_i( k)⟩=0,     k∈,  i=0,1,2,...,d.

The above characterization can also be understood that on BISs, the resonant spin-reversing transitions γ⃗→-γ⃗ are induced by the SO field h_ so( k) due the vanishing gap (for band-crossing surfaces). Thus the time-averaged spin-polarizations vanish at BISs. In contrast, away from the BISs, the nonzero h_0( k) severs as a detuning for the spin-reversing transitions, so the time-averaged spin-polarizations do not vanish.

The vanishing ⟨γ_i⟩ on BIS implies that the topological number of the dD gapped phase cannot be obtained via the time-averaged spin-texture ⟨γ⃗( k)⟩. Interestingly, more nontrivial features of the quench dynamics are captured by the variation of the time-averaged spin polarizations across the BISs. Note that the component h_0<0 in the vector volume V_ BIS and h_0>0 in the volume V̅_ BIS. When passing through a BIS, the relative angle ϕ between the initial spin-polarization and the vector h⃗( k) varies from ϕ<π/2 in V_ BIS, with a nonzero dynamical averaging ⟨γ⃗⟩ pointing to h⃗( k)-direction, to ϕ>π/2 in V̅_ BIS, with nonzero ⟨γ⃗⟩ pointing oppositely to h⃗( k)-direction (Fig. 1c). This follows that the variation of ⟨γ⃗⟩ across the BIS follows the direction of the spin-orbit field h_ so. To quantify this picture, we define a new dynamical spin-texture field g( k), whose components are given by g_i( k)≡-1/ N_ k∂_k_⊥⟨γ_i⟩, with N_ k being the normalization factor. Here k_⊥ denotes the momentum perpendicular to BIS and points from V_ BIS to V̅_ BIS. It can be shown from Eq. (<ref>) that on the BISs

    g_i( k)|_ k∈=ĥ_ so,i( k).

With this identity we reach immediately that the topological invariant defined on the BISs can be recast into a dynamical form, namely

    w_d-1=∑_jΓ(d/2)/2π^d/21/(d-1)!∫_BIS_j g( k)[ g( k)]^d-1.

The right hand side of the above formula describes the mapping from BISs to (d-1)D spherical surface which the emergent dynamical field g( k) belongs to, rendering the dynamical classification of the dD gapped phase. The dynamical topological invariant is intuitively described by the coverage of the dynamical spin-texture g( k) over the full (d-1)D spherical surface.
The results in Eqs. (<ref>)-(<ref>) manifest a highly nontrivial dynamical bulk-surface correspondence, which can be directly measured from the dynamical spin-polarization patterns emerging on the BISs. Since the quantum spin dynamics is resonant and nontrivial only on BISs, the dynamical bulk-surface correspondence can be well resolved in experiments. This enables us to propose high-precision dynamical schemes to detect the bulk topology, as studied below with experimentally feasible models.

We note that the dynamical classification theory is essentially based on the mapping of the dD bulk topology to (d-1)D BISs. This feature is valid for generic multiband systems, in which case the BISs are defined among every n_d bands in the BZ. The topology only relies on the (d-1)D invariants `locally' defined on the BISs, which can be measured based on a sequence of quenches by initializing the system in different bands (Appendix). This dimension reduction is the essential difference from the previous dynamical schemes which are applicable to two-band models <cit.>.

Application to quenching experiments.–We first consider the 1D topological phases of AIII class obtained by the Hamiltonian ℋ(k_x)=h⃗(k_x)·σ⃗, where h_x=t_sosin k_x≡ h_ so, h_y=0, and h_z=m_z-t_0cos k_x≡ h_0, with t_0 and t_ so denoting the nearest-neighbor spin-conserved and spin-flipped hopping coefficients. The model was proposed in Ref. <cit.> and realized recently in experiment for ^173Yb fermions <cit.>. The effective magnetization m_z can be precisely tuned by a bias magnetic field. The quench study is performed by tuning the Zeeman term from
m_z≫ t_0, which corresponds to the deep trivial
regime with the system initialized in spin-down, to |m_z|<t_0, which lies in topological regime.
The quantum dynamics is then governed by the unitary evolution operator (see Methods).
The numerical results are shown in Fig. 1d-e, with the parameters taken as t_so=0.2t_0. The resonant spin-flip transitions driven by SO field emerge at two momentum points k_ BIS=±cos^-1(m_z/t_0) with vanishing time-averaged spin polarizations ⟨σ⃗⟩, which manifests the BIS (band inversion points). Away from the BIS the spin-flip oscillation is well suppressed by the energy gap. The opposite magnitudes of the directional derivative ∂_k_⊥⟨σ_x⟩ at two band inversion points show the nontrivial topology of the post-quench regime, and gives the winding number +1 for |m_z|<t_0 (Fig. 1e).





In 2D case we consider the quantum anomalous Hall model <cit.> realized in recent experiments <cit.> ℋ(𝐤)=h⃗(𝐤)·σ⃗,
where the vector field reads h⃗(𝐤)=(m_x+t_so^xsin k_x,m_y+
t_so^ysin k_y,m_z-t_0cos k_x-t_0cos k_y).
The flexibility in the decomposition of the vector h⃗(𝐤) allows us to consider different quench ways, e.g., by quenching h_z or h_y axis, which can benefit real experimental measurement.
First, we take h_0≡ h_z,  h_ so≡(h_x,h_y), and the quench is performed by varying m_z suddenly from m_z≫ t_0 to 0<|m_z|<2t_0, while setting m_x=m_y=0. Similar to 1D model, the resonant transition between the spin-down and spin-up states only occurs at the BIS (band inversion ring). From the numerical results given in Fig. 2a, with m_z=t_0 and t_ so^x,y=0.2t_0, the vanishing time-averaged spin-polarizations ⟨σ_i⟩=0 are obtained on the BIS. Further, the emergent dynamical field g( k)=-∂_k_⊥⟨σ⃗⟩ exhibits a nonzero 1D winding on the BIS, which characterizes a topological charge located at Γ point [with k=(0,0)] and gives the Chern number C_1=-1 (Fig. 2b). Similarly, if -2t_0<m_0<0 for the post-quench regime, an opposite dynamical winding exhibits on the BIS which circles M point [with k=(π,π)], and gives C_1=+1 (see Appendix). Secondly, we take h_0≡ h_y and quench the h_y axis by ramping m_y quickly from m_y≫ t_0 to m_y=0, while setting m_x=0 and m_z=t_0. The time-averaged dynamical spin textures exhibit line-shape structures (Fig. 2c), which correspond to the two open BISs at k_y=0 and -π, respectively. Despite a trivial pattern along k_x with k_y=-π, the winding of the dynamical spin-texture g( k) (in x-z plane) along k_x with k_y=0 measures the Chern number C_1=-1 (Fig. 2d). Applying the dynamical classification to measuring topological states with high Chern numbers is straightforward and is found in Appendix.




We further consider the application to 3D topological phases, whose hosting Hamiltonian reads ℋ(𝐤)=h⃗(𝐤)·γ⃗,
with h_0(𝐤)=m_z-t_0∑_icos k_i and
h_i=t_sosin k_i (i=x,y,z). Here we take that
γ_0=σ_z⊗τ_x, γ_1=σ_x⊗𝕀,
γ_2=σ_y⊗𝕀 and γ_3=σ_z⊗τ_z, where the Pauli matrices
σ_i and τ_i may refer to the real- and pseudo-spin degrees of freedom (e.g., sublattices or orbitals). The Hamiltonian has a chiral symmetry defined by σ_z⊗τ_y, so it belongs to AIII class and is classified by 3D winding numbers in equilibrium theory <cit.>. The trivial phase corresponds
to |m_z|>3t_0, while the topological phases include three regions:
(I) t_0<m_z<3t_0 with winding number ν_3=-1; (II) -t_0<m_z<t_0
with ν_3=2; and (III) -3t_0<m_z<-t_0 with ν_3=-1.
We consider the quench from m_z≫ t_0 to say phase I (m_z=1.3t_0 and t_so=0.2t_0), with the numerical results being partially presented in Fig. 3 (more is shown in Appendix).
The BIS of the post-quench band is given in Fig. 3a, with the vector arrows denoting the dynamical field g( k). In a practical measurement, by resolving the time-averaged spin polarizations ⟨γ⃗⟩ on two closed surfaces sightly inside (with h_0=-0.1t_0) and outside (with h_0=0.1t_0) the BIS, one can qualitatively determine the vector field g( k) (without affecting measurement of topology) from the subtraction of ⟨γ⃗⟩ between the two surfaces. A few examples for ⟨γ_3⟩ are given in Fig. 3b-d, similar for ⟨γ_1,2⟩. From the dynamical field g( k) on the BIS, as illustrated with arrows in Fig. 3a, which gives a nonzero Chern number C_1=-1 on the 2D BIS, corresponding to the 3D winding number ν_3=-1 of the post-quench 3D system.




While the above study is based on unitary quench dynamics, the dynamical classification can be generalized to the dissipative quantum dynamics at finite temperature. In the presence of dissipation, we expect that quantum dynamics should be nearly unitary within a short-time evolution after quench, while the dissipation effects shall dominate in the long-term evolution. Taking the 2D model as an example,
we describe the dynamics by the Lindblad master equation <cit.>

    / tρ_𝐤=-[ℋ,ρ_𝐤]+η(σ̃_-ρ_𝐤σ̃_+
    -1/2{σ̃_+σ̃_-,ρ_𝐤})

with initial state being ρ_𝐤(0)=f(E_+,T)|+,𝐤⟩⟨ +,𝐤|
+f(E_-,T)|-,𝐤⟩⟨ -,𝐤|.
Here the distribution function f(E,T)=[e^(E-μ)/k_BT±1]^-1 if the system is simulated with fermions (for +) and bosons (for -), |±,𝐤⟩
are the eigenstates of the pre-quench Hamiltonian,
the Pauli matrices σ̃_±≡(σ̃_x±σ̃_y)/2
are defined in the eigenbasis of the post-quench Hamiltonian
ℋ, and η is the decay
rate. A numerical study with typical parameters is shown in Fig. 4. For a short-term evolution, we find that the dissipative effect is negligible  (nearly unitary), as demonstrated in Fig. 4a (gray area) and b. In this regime, the previous dynamical classification is well applicable to measure the locations of BISs and the bulk topology (Fig. 4b).
Further, after a long-term evolution, the dissipation
drives the occupation of Bloch states to approach equilibrium results. In this case we find that the spin-polarization ⟨σ⃗( k)⟩ gradually points oppositely to the h_ so( k)-direction on BISs of the post-quench system (Fig. 4a). Then the Chern number can be directly read out by the long-term spin texture ⟨σ⃗( k)⟩, without the further adoption of g( k) for the classification (Fig. 4a,c).

Discussion and outlook.–The dynamical classification exhibits essential advantages in experimental investigation of topological quantum states. In particular, applying the classification theory proposed in this work to measuring 2D Chern insulators has been achieved in a very recent experiment <cit.>, which confirms that observing topological phases based on the current dynamical classification is of much higher precision, compared with that based on equilibrium classifications. The advantages are rooted in two essential aspects. First, the bulk-surface duality uncovered here maps the classification of bulk topology to lower-dimensional invariants on BISs, and simplifies the topological characterization. Secondly, the spin dynamics is resonant and nontrivial only on the BISs, so the dynamical bulk-surface correspondence can be easily resolved in the quench study. Finally, the dynamical scheme by nature is robust against the non-ideal conditions in both the state initialization and measurement. In particular, the system starts at deep trivial regime, with the initial phase being dependent on only a single parameter, Zeeman term, immune to non-ideal conditions. Quenching it to topological regime induces quantum dynamics which exhibits resonant oscillations on the BISs. The short-term unitary quench dynamics characterizing the topology of the post-quench system is not affected by detrimental effects like the thermal effects, leading to high-precision measurement of full topological phase diagrams, as confirmed in experiment <cit.>. In comparison, measuring the topology of static phases necessitates the systems to be carefully initialized, with the quantum states being properly occupied. The preparation is intrinsically sensitive to non-ideal conditions including thermal effects, especially unsatisfactory in the regime close to phase boundaries <cit.>. Our dynamical classification provides new approaches with high feasibility in exploring topological quantum physics.

Moreover, the dynamical classification does not require the initial trivial phase to be fully polarized along certain direction. In the generic case with a trivial initial phase which has trivial but momentum-dependent (pseudo)spin texture, one can rotate the initial (pseudo)spin of each Bloch state to the direction of h_0 axis by a momentum-dependent unitary transformation. Accordingly, under the transformation the new (pseudo)spin bases become momentum-dependent but still topologically trivial. Then, if expressed in the new (pseudo)spin bases, all the results including the characterization of BISs by vanishing time-averaged spin polarization and dynamical bulk-surface correspondence take the same forms as in the case with a fully polarized initial phase. However, taking the fully polarized (deep trivial) phase as the initial state is most convenient for real experiments.

Finally, this work provides important insights into developing non-equilibrium classification for the complete classes of equilibrium topological quantum states, as conceptually different from the known various equilibrium classification theories. The dynamical classification theory established here covers different categories of the topological phases characterized by integer invariants, including the A, AIII, BDI, C, CII, and D classes in the AZ ten-fold ways <cit.>. It is promising and of great interests to generalize the current classification theory to all ten-fold categories <cit.>, crystalline topological states, including the insulators and superconductors protected by space groups <cit.>, and further to the correlated topological states of fermions and bosons <cit.>. We note that for the correlated topological phases whose topology can be characterized by quasiparticles or mean-field picture, the concept of BISs can still be well defined and the present dynamical classification theory should be applicable. In turn, the correlation effects can have non-trivial effects on the dynamical bulk-surface correspondence and deserve in-depth studies.. We therefore expect that the present work can open a broad new direction to classify and detect topological quantum states by non-equilibrium quantum dynamics.


This work was supported by the National Key R&D Program of China (2016YFA0301604), National Nature Science Foundation of China (under grants No. 11574008 and No. 11761161003), and the Thousand-Young-Talent Program of China.












§ APPENDIX A: CHARACTERIZING DD TOPOLOGICAL PHASES BY (D-1)D INVARIANTS


  


In this section, we present the rigorous proof of the bulk-surface duality that dD topological phases can be characterized by
the (d-1)D topological invariant defined on the BISs, including the generalization to multiband systems.




 §.§ 1. Bulk-surface duality: the (d-1)D topological invariant on the BISs





For a generic dD gapped topological phases described by the basic Hamiltonian (<ref>),
there can be multiple BISs defined by BIS≡{𝐤|h_0( k)=0} in the BZ, being closed to connected BZ boundary (Fig. <ref>a-c). For convenience, here we denote
by 𝒱_ BIS the region enclosed by BISs with h_0(𝐤)<0,

and decompose the momentum 𝐤 into (k_⊥,k_∥),
with k_⊥ being perpendicular to the BIS and pointing to the region 𝒱̅_BIS where h_0(𝐤)>0.

Note that the SO field opens the gap at the BISs. We shall show that the SO field on the BISs determines the topology
uniquely. To prove it, we take advantage of the fact that the topology remains
unchanged under any continuous deformation that does not close the gap <cit.>. Therefore,
we can keep the SO field only in the region near the BISs but turn it off
elsewhere, which results in the following deformed Hamiltonian

    ℋ(𝐤)→h⃗'(𝐤)·γ⃗=sinΘ_𝐤γ_0+cosΘ_𝐤∑_i=1^dĥ_so,i(𝐤)γ_i,

where ĥ_so(𝐤)≡h_so(𝐤)/|h_so(𝐤)| denotes the unit SO field. Here Θ_𝐤∈[-π/2,π/2] is a smooth and
monotonic function of h_0(𝐤), which approaches the step function and maps
h_0(𝐤) to ±π/2 if h_0(𝐤)≷0
or to 0 if h_0(𝐤)=0, see Fig. <ref>d.
For example, it can be chosen as the error function πerf(ξ h_0(𝐤))/2
with ξ→∞. On the other hand, the
SO field with pre-factor cosΘ_𝐤 is truncated to the region near the BISs. It is noteworthy
that this deformation keeps the information of BISs including the SO field on
it, and the deformed Hamiltonian has been normalized.

The bulk-surface duality can be explicitly proven with the deformed Hamiltonian (<ref>).
The proof is given as follows. Generally, the topological phases in d=2n-1 dimensions are characterized
by the (2n-1)D winding number <cit.>

    N_2n-1=(-1)^n-1(n-1)!/2(2π)^n(2n-1)!∫_ BZ Tr[γℋ(ℋ)^2n-1],

where γ=()^n∏^2n-1_i=0γ_i is the chiral symmetry for the Hamiltonian, satisfying γ^†=γ, γ^2=1 and {γ,γ_i}=0.
By using the deformed Hamiltonian (<ref>), we obtain


    Tr[γℋ(ℋ)^2n-1]=(2n-1)!(-2)^n∑_i=0^2n-1ϵ^i01⋯i̅⋯(2n-1)h'_i h'_0∧ h'_1∧⋯∧ h'_i̅∧⋯∧ h'_2n-1,


where `d' denotes the exterior derivative, i̅ means that the index does not take the value i. Here the trace property of Clifford matrices Tr[γγ_iγ_j⋯γ_sγ_t]=(-2)^nϵ^ij⋯ st for dimensions d=2n-1 (see Supplementary Material), contributes the factor (-2)^n, while the product of the Levi-Civita tensor and the wedge will lead to the term (2n-1)!ϵ^i01⋯i̅⋯(2n-1)=(2n-1)!(-1)^i. Similarly, in d=2n dimensions, via the the trace property of Clifford matrices Tr[γ_iγ_j⋯γ_sγ_t]=(-2)^nϵ^ij⋯ st (see Supplementary Material), the topological phases are characterized by the dD (or the n-th) Chern number <cit.>

    Ch_n=-1/2^2n+11/n!(/2π)^n∫_ BZ Tr[ℋ(ℋ)^2n],

with


    Tr[ℋ(ℋ)^2n]=(2n)!(-2)^n∑_i=0^2nϵ^i01⋯i̅⋯(2n)h'_i h'_0∧ h'_1∧⋯∧ h'_i̅∧⋯∧ h'_2n.

From the above result, one can observe that the winding number and Chern number have the same structure and are proportional to the terms ∑^d_i=0(-1)^ih'_i⋀̅^d_l=0 h'_l, where the wedge product ⋀̅ indicates the index l≠ i. It is therefore convenient to consider these terms at first, which can be readily calculated with the deformed Hamiltonian (<ref>):










    ∑^d_i=0(-1)^ih'_i⋀̅^d_l=0 h'_l=sinΘ_𝐤cos^dΘ_𝐤⋀_l=1^dĥ_so,l+cos^d-1Θ_𝐤Θ_𝐤∧∑_i=1^d(-1)^iĥ_so,i⋀̅_l=1^dĥ_so,l.


To further simplify the above formula, one can choose a special configuration for the Θ_ k function, while the final result is independent of the choice of Θ_ k. For this we let the SO field (∝cosΘ_𝐤) for the deformed Hamiltonian be truncated to the BISs (sinΘ_𝐤=0) by taking Θ_𝐤 to
approach the step function. Then the contribution to the integral from the first term on the r.h.s. of Eq. (<ref>) vanishes. On the other hand, since Θ_𝐤 tends to be a delta function, the whole integral can be reduced to that over the region ℛ with non-zero Θ_𝐤 (see Fig. <ref>c-d). By decomposing the measure ^d𝐤 into k_⊥^d-1k_∥ (see Fig. <ref>c), the integral for Eq. (<ref>) over the whole BZ can further reduce to the integral over the BISs, as shown by


    ∫_BZ∑^d_i=0(-1)^ih'_i⋀̅^d_l=0 h'_l
    
    =    ∫_ℛcos^d-1Θ_𝐤Θ_𝐤∧∑_i=1^d(-1)^iĥ_so,i⋀̅_l=1^dĥ_so,l
    
    =    ∫_-π/2^π/2Θ_k_⊥cos^d-1Θ_k_⊥
           ×∑_j∫_ BIS_j∑_i=1^d(-1)^iĥ_so,i⋀̅_l=1^dĥ_so,l.


















Finally, by a straightforward calculation of the integral ∫^π/2_-π/2 x cos^d-1x, both the winding number and Chern number can be expressed in a unified form as the winding of the SO field ĥ_so along the BISs:

    w_d-1=∑_jΓ(d/2)/2π^d/21/(d-1)!∫_ BIS_jĥ_so(ĥ_so)^d-1,

where ĥ_so(ĥ_so)^d-1≡ϵ^i_1i_2⋯ i_dĥ_so,i_1ĥ_so,i_2∧⋯∧ĥ_so,i_d is defined on the (d-1)D BISs, with the indices i_1,…,d∈{1,…,d}, and Γ(x) is the Gamma function with Γ(x+1)=xΓ(x),
Γ(1)=1 and Γ(1/2)=√(π).
Note that the prefactor 2π^d/2/Γ(d/2) is actually the area of a (d-1)D unit sphere.
Analogous to the “Gauss's law”, the invariant w_d-1 simply counts the total topological charges of the SO field located at momenta with h_so( k)=0 and enclosed by the BISs. Finally we reach

    ν_d=w_d-1=∑_i∈𝒱_BIS𝒞_i,

where

    𝒞_i≡Γ(d/2)/2π^d/21/(d-1)!∫_BIS_iĥ_ so(ĥ_ so)^d-1

represents the total topological charge in the patch of V_ BIS enclosed by the i-th BIS. The topological phase transition occurs when a charge passes through a BIS, since the gap must be closed and reopen as the topological charge located at h_so( k)=0 passes through the BISs with  h_0( k)=0. Accordingly, the topology characterized by the SO field on BISs is changed.



To gain intuition, we take the Chern insulator as an example, with Hamiltonian given by ℋ(𝐤)=h⃗(𝐤)·σ⃗ with h⃗(𝐤)=(t_sosin k_x,t_sosin k_y, m_z-t_0cos k_x-t_0cos k_y). According to the bulk-surface duality, the Chern number reduces to the form Ch_1=∑_j1/2πϵ^mn∮_BIS_jℓ̂·ĥ_so,m∇ĥ_so,n,
counting the topological charges of the SO field h_so within V_ BIS <cit.>.
First, we take that h_0≡ h_z and h_so≡ (h_y,h_x). The topological charges and BISs are shown in Fig. <ref>a, which explicitly determine Chern number for different phases governed by m_z: The Chern number Ch_1=-1 for 0<m_z<2; Ch_1=+1 for -2<m_z<0; Ch_1=0 for otherwise.
Further, if we choose h_0≡ h_x such that 𝒱_BIS
with h_0(𝐤)<0 corresponds to the region with -π<k_x<0, and then
the SO field is h_so≡(h_z,h_y).
Two topological charges with opposite signs, located at 𝐤=(±arccos[(m_z+1)/t_0],π) for -2<m_z<0
or at k=(±arccos[(m_z-1)/t_0],0) for 0<m_z<2, are created from vacuum,
moves and finally annihilates, as m_z varies from -∞ to
+∞ [see Fig. <ref>b], which characterizes the same Chern number for each topological regime as in Fig. <ref>a.



 §.§ 2. Generalization to generic multiband systems


We show here that the bulk-surface duality is valid for the generic systems with arbitrary multiple
bands and multiple band crossings in the BZ, as illustrated Fig. <ref>c. Note that for a generic dD topological phase, the gap must be opened through each group of n_d bands (n_d=2^d/2 or 2^(d+1)/2) with hybridization at BISs of these bands.
Therefore, it is more convenient, by transforming the system into a new trivial bases, which block diagonalize the Hamiltonian ℋ(𝐤)=ℋ_1(𝐤)⊕ℋ_2(𝐤)⊕⋯ at each momentum 𝐤. Only those blocks involving n_d bands with band crossings have nontrivial contribution to the topology. Note that this transformation does not change the topology of the Block Hamiltonians due to the triviality of the new bases.
For convenience, we denote each band crossing for a group of n_d bands as a BIS. Without changing the topology of the system, we assume that each BIS is separated from others in the BZ (If there are more than one BIS occurring at the same k point, one can continuously deform the band dispersion and separate such BISs without affecting the bulk topology). In this case, the Hamiltonian around each BIS can be decomposed to the contributions from three parts: the n_d bands with hybridization on the BIS, the bands with energy above the BIS, and the bands below the BIS. Thus the Hamiltonian can be locally written as ℋ(𝐤)=𝒫_upper(𝐤)⊕ℋ_BIS(𝐤) ⊕𝒫_lower(𝐤).
Here 𝒫_upper/lower are the projection
operators for the bands with higher (lower) energy than the part with
band crossing ℋ_BIS. For the higher or lower bands 𝒫_upper/lower, as there is no band crossing, the corresponding spin textures are trivial and no singularities exist for
the corresponding eigenstates. One can safely deform the bands of 𝒫_upper/lower
into constants in dispersion. Similarly, away from the BISs, no band crossing occurs, so all the bands can be smoothly connected in those areas of the BZ, and also be further deformed into constants in dispersion. It is straightforward to know that the flattened pathes in BZ have no contribution to the bulk topology, which is then solely determined from the integral over the BISs. Further, to compute the contribution from a specific BIS, we can single out the n_d bands which are characterized by the Hamiltonian ℋ_BIS(𝐤) on the BIS, and flatten these bands in all the remaining k space of BZ away from such BIS. This can be safely done since the band crossings occurring at other momenta of BZ are unrelated to the contribution to bulk topology from the BIS under consideration (thus those band crossings can be removed for this computation). In this way, the contribution of the BIS is effectively determined by a basic Hamiltonian with the band crossing at BIS described by ℋ_BIS(𝐤), precisely given by the (d-1)D topological number defined on the BIS. The total topological invariant of the dD phase is then a summation of contributions from all BISs, formed by all bands of the BZ.
This completes the proof.

The above proof implies that the (d-1)D topological number shown in Eq. (<ref>) is also valid for the case with multiple bands. However, in this generic case the SO field h_ so ( k) is defined locally on each BIS, characterizing only the hybridization among the n_d bands of such BIS.



§ APPENDIX B: DYNAMICAL BULK-SURFACE CORRESPONDENCE


  

In this section, we present the details of the dynamical bulk-surface correspondence.
Two basic issues need to be addressed: how to identify the BISs and how to extract the dynamical topological invariant.
The study is first based on the basic Hamiltonian (<ref>), and then generalized to generic multiband systems.



 §.§ 1. Characterizing the BISs by quench dynamics


As discussed in the main text, the
system is initially prepared in the ground state of the pre-quench Hamiltonian
ℋ(m_0|_t<0,𝐤) in the deep trivial regime, with density matrix ρ_0, where the spins are
fully polarized in the opposite direction of the axes h_0.

After a sudden quench to a topological
non-trivial Hamiltonian ℋ(m_0|_t>0,𝐤), each momentum-linked spin
will rotate around the post-quench vector field h⃗(𝐤), with the quantum dynamics being governed by the unitary evolution operator U(t)=exp(-ℋt) for the post-quench system, see Fig. <ref>c.







We consider the time-averaged spin polarizations, which can be measured in experiments:

    ⟨γ_i⟩=lim_T→∞1/T∫_0^T t Tr[ρ_0 e^ℋtγ_ie^-ℋt].

For the Hamiltonian (<ref>), one can show that e^ℋt=cos(Et)+sin(Et)ℋ/E with E( k)=√(∑_i=0^dh_i^2) the energy. Thus the time-averaged spin polarizations read

    ⟨γ_i( k)⟩=h_iTr[ρ_0ℋ]/E^2=-h_i( k)h_0( k)/E^2( k).

On the BISs with h_0(𝐤)=0 of the post-quench Hamitlonian, one can find that all the time-averaged spin polarizations ⟨γ_i( k)⟩ (for i=0,1,...,d) vanish on the BISs.
Hence, the BISs can be dynamically determined by the surface with vanishing time-averaged spin polarizations:

    BIS={𝐤|⟨γ⃗(𝐤)⟩=0}.





 §.§ 2. Dynamical topological invariants and the dynamical bulk-surface correspondence


Time-averaged spin polarizations ⟨γ_i⟩ vanish on the BISs, but are nonzero inside and outside the BISs. The variation of ⟨γ_i⟩ across the BISs can characterize the topology. For this a new dynamical spin-texture field is defined in the main text as the directional derivative of time-averaged spin textures ⟨γ⃗⟩ across the BISs:

    g(𝐤)   ≡-1/ N_ k∂_k_⊥⟨γ⃗⟩=-lim_k_⊥→ 0^+1/ N_ kΔ⟨γ⃗⟩/2k_⊥
       =-lim_k_⊥→ 0^+1/ N_ k⟨γ⃗⟩|_h_0(k_⊥)>0-⟨γ⃗⟩|_h_0(-k_⊥)<0/2k_⊥,

where N_ k is the normalization factor, and the origin of k_⊥ has been shifted locally to the BISs.
Note that for the pre-quench Hamiltonian with h_0(𝐤)≪ 0, the time averaged spin polarizations ⟨γ_i⟩=+h_ih_0/E^2. Accordingly, the dynamical spin-texture field should be defined by g(𝐤)≃∂_k_⊥⟨γ⃗⟩, which takes an additional minus sign.

Near the BISs, the dispersion h_0(𝐤) is monotonous and is determined by k_⊥ uniquely. One can identify h_0 from k_⊥, that is, h_0≅ k_⊥ up to some local scale transformation that does not affect the topology. The difference of ⟨γ_i⟩ between the two sides of BISs reads


    Δ⟨γ_i⟩|_k_⊥→0   =    - [h_i(k_⊥,k_∥)k_⊥/k_⊥^2+∑_j=1^dh_j^2(k_⊥,k_∥)-
          -h_i(-k_⊥,k_∥)(-k_⊥)/k_⊥^2+∑_j=1^dh_j^2(-k_⊥,k_∥)]
       =   - 2h_i(0,k_∥)k_⊥/∑_j=1^dh_j^2(0,k_∥)+𝒪(k_⊥^2).


We obtain the emergent dynamical spin-texture field by

    g_i(𝐤)|_ k∈ BISs=h_i(0,k_∥)/∑_j=1^dh_j^2(0,k_∥)=ĥ_so,i.

With this result we reach that the bulk topology can be characterized by the dynamical invariant defined as the winding of the new dynamical field g(𝐤) on the BISs:

    w_d-1≡∑_jΓ(d/2)/2π^d/21/(d-1)!∫_ BIS_j g(𝐤)[ g(𝐤)]^d-1.

This is the dynamical bulk-surface correspondence.



 §.§ 3. Quench dynamics in multi-band systems


Now we show the dynamical classification for multi-band systems. To better understand our results, we consider first the weak spin-orbit coupling regime, where the band dispersions are dominated by the bare energies ϵ_j( k) of each bare orbital λ_j which define the BISs, and the spin-orbit couplings are weak and induce small gap on the BISs. Then we consider the generic case, including the strong spin-orbit coupling regime, namely, the spin-orbit couplings can even dominate the dispersions of the topological bands. For convenience, we consider here the 1D or 2D topological phases, which are most relevant for the current cold atom experiments, and the proof can be directly applied to the higher dimensional phases.

For the weak coupling regime, the spin-orbit couplings are weak, so that they only have local effects in the
momentum space and cannot change the band structure determined by the bare energies of orbitals dramatically. Note that for a 1D/2D system, a key feature is that the topological gap must be opened between every two bands. Therefore, all one needs to determine is the information for all these band crossings determined by these bare orbitals, namely, the BISs defined between every two bands. Then the detection scheme for the minimal Hamiltonian can be generalized directly. In
particular, we initialize the system in the deep trivial regime, with all subbands for orbitals being far separated from each other (similar to the large Zeeman splitting in the minimal model) and prepare the system in a single band (denoted as λ). Then we pick up this band and another one λ_1, and quench the effective Zeeman splitting m_1 between the two bands (the energy splitting between them) to desired post-quench regime, but the parameters for all the remaining bands keep unchanged. With this process we can dynamically determine the BISs defined by λ and λ_1 bands, and the topological invariants defined on the BISs, by measuring the Pauli operators defined in these two orbitals. Then, replacing λ_1 with another band λ_2 and through a similar process we can detect the topology on the BISs defined by λ and λ_2 bands. With λ_j running over all the bands, we obtain the topology of λ band by summing over the invariants detected in these processes. Repeating this detection for all λ bands we obtain the topology for all subbands of the whole system. It is trivial to know that when the band number is two, a single quench determines the topology, reducing to the minimal model regime.

Now we show the detection scheme for the generic case without assuming the spin-orbit coupling to be weak. The key difference between the weak and strong spin-orbit coupling regimes is that, in the strong coupling regime the band crossing between two subbands can be affected by the couplings between such two bands and the remaining bands. Thus picking up only every two subbands for quench may not extract the true topology. Nevertheless our detection scheme can be still applied. Again, the key observation is that the topological gap has to be opened between every two bands (for 1D/2D system). Therefore, at every momentum k, by transforming the system into a new trivial bases, we can always block diagonalize the Hamiltonian ℋ( k)=ℋ_1(𝐤)⊕ℋ_2(𝐤)⊕⋯, with each block ℋ_j being at most two by two matrix.  More generically, this block diagonal picture is necessary only in the vicinity of each BIS.
Note that all the band crossing information can be kept in each block diagonalized Hamiltonian ℋ_j, and only these blocks with band crossings have nontrivial contribution to the topology, with the topology determined by the corresponding BISs defined by h_0^(j) axis in ℋ_j. All we need to do is detect the topological invariants on these BISs. However, in the present block diagonal picture, the (pseudo)spin-1/2 bases for each block Hamiltonian ℋ_j is in general momentum dependent, rather than being constant, while such bases do not have nontrivial topology. This is the key difference from the minimal model and weak coupling regime. The dynamical detection of the topology can be performed in the following way. We initialize the system in the deep trivial regime with the energies of the bare orbitals being far separated from each other. The bare orbitals are again denoted as λ_j, and note that they are momentum-independent. The (pseudo)spin bases of each block Hamiltonian ℋ_j, denoted by ↑^(j)_ k and ↓^(j)_ k, are obtained by momentum-dependent unitary transformations on the bare orbitals. The system is prepared in certain orbital λ, and is then quenched to the final total Hamiltonian ℋ. For each block ℋ_j, the initial state is simply obtained by projecting λ-state into the (pseudo)spin bases of ℋ_j, i.e., |0_ k⟩=(c_↑(𝐤),c_↓(𝐤))^𝖳, which is momentum-dependent. This is equivalent to the case that for the minimal model, the initial trivial state is not fully polarized, but momentum-dependent. Accordingly, in quench dynamics, we can measure the momentum-dependent observable σ^(j)_z; k=|0_ k⟩⟨ 0_ k|-|0'_ k⟩⟨ 0'_ k|, where |0'_ k⟩ is the state perpendicular to |0_ k⟩, and of which the time-averaged expectation value reads

    ⟨σ^(j)_z; k⟩    =lim_T→∞1/T∫^T_0 t Tr(|0_ k⟩⟨ 0_ k|e^ℋ_j( k)tσ^(j)_z; ke^-ℋ_j( k)t) 
        = lim_T→∞1/T∫^T_0 t Tr(|0⟩⟨ 0|e^ℋ'_j( k)tσ^(j)_ze^-ℋ'_j( k)t),

where the second line is obtained by a unitary transformation U(𝐤)=(|0_ k⟩,|0'_ k⟩), and the Hamiltonian ℋ'_j(𝐤)=U( k)ℋ_j(𝐤)U^†( k) is topologically equivalent to ℋ_j since |0_ k⟩ and |0'_ k⟩ are trivial. Here |0⟩=(1,0)^𝖳 and σ^(j)_z is the usual Pauli matrix. The above equality means that one has chosen a momentum-dependent σ^(j)_z; k axis to define the BISs for ℋ_j block, which coincide with the BISs for the Hamiltonian ℋ'_j with the same topology as ℋ_j. Therefore, this momentum-dependent h^(j)_0 axis does not affect the existence of BISs, nor the topology defined on such BISs.
Then, as long as the projection of λ-state to the (pseudo)spin bases of ℋ_j is nonzero, from the (pseudo)spin dynamics one can extract the dynamical topological invariant for the block Hamiltonian ℋ_j. Further, by initializing the system in different bare orbital λ-states and performing quench until the projection of initial states can cover (pseudo)spin bases of all Hamiltonian blocks, one can obtain the topological invariants for all ℋ_j blocks, which characterize the topology of all bands of the total Hamiltonian.





We note that in the real cold atom experiment, one can always simplify the measurement by well controlling the parameters, e.g. staring from the weak spin-orbit regime. In this way, one can get enough information of the bands and topology before entering the strong coupling regime, and simplify the detection in generic case. Therefore, our dynamical classification theory is universal for the case of multiband systems.



§ APPENDIX C: NUMERICAL RESULTS


In this section, we provide more numerical calculations related to the main text for the quench study.



 §.§ 1. Topological phases with high invariant


The dynamical classification theory is clearly applicable to topological phases with high invariants. To illustrate the validity, we consider two different 2D models H_j( k)=h⃗_j( k)·σ⃗ (j=1,2), which are variations of the 2D QAH model studied in the main text.

We consider the first model H_1 with the vector field h⃗_1( k)=(t_ so
sin 2k_x,t_ sosin2k_y,m_z-t_0cos k_x-t_0cos k_y). While the trivial phase corresponds to |m_z|>2t_0, the topological phases are distinguished as: (i) t_0<m_z<2t_0
with the Chern number C_1=-1; (ii) 0<m_z<t_0 with C_1=3; (iii)  -t_0<m_z<0 with C_1=-3; (iv) -2t_0<m_z<-t_0 with C_1=1.
When t_0≫ t_ so, this model has almost the same band structure as the 2D model in the main text.
It is then supposed that
the results of quenching to the phases with C_1=±1 should be similar to those shown in Fig. 2a. Here we only provide the numerical
results of quenching to the phase with C_1=3 (see Fig. <ref>a-b), where the time-averaged spin textures also exhibit ring-shape structure, but
the corresponding spin-texture field winds around three times, charactering the high Chern number.

We then consider the other one H_2 with h⃗_2( k)=(t_ sosin 2k_x,t_ sosin2k_y,m_z-t_0cos 2k_x-t_0cos 2k_y). Apart from the trivial phase with |m_z|>2t_0, this model has two topological phases:
(i) 0<m_z<2t_0 with the Chern number C_1=-4; (ii) -2t_0<m_z<0 with C_1=4. After a sudden quench from m_z≫ 2t_0 to |m_z|<2t_0, the time-averaged spin textures have four rings (see Fig. <ref>c-d), and the dynamical spin-texture field on each ring winds around once. Thus four rings in the BZ characterize the Chern number C_1=±4. More results, including the case with a trivial phase but having multiple BISs, are discussed in Supplementary Material <cit.>.




 §.§ 2. Characterizing 3D topological phases


As discussed in the main text, the 3D model H=h⃗( k)·γ⃗, which is a direct generalization of the 2D model (see the main text), has three topological phases. In Fig. <ref>, we have shown the BIS of the topological phase with 3D winding number ν_3=-1 and the time-averaged spin texture
⟨γ_3⟩. Here we provide more results and demonstrate that the dynamical spin-texture field g(𝐤) on the BIS is qualitatively the same as that determined from the difference of ⟨γ_j⟩ (j=1,2,3)
on the two equal-energy surfaces slightly outside (h_0=0.1t_0) and inside (h_0=-0.1t_0) the BIS (see Fig. <ref>).



99
Klitzing1980 K. v. Klitzing, G. Dorda, M. Pepper, New Method for High-Accuracy Determination of the Fine-Structure Constant Based on Quantized Hall Resistance. Phys. Rev. Lett. 45, 494 (1980).
Tsui1982 D. C. Tsui, H. L. Stormer, A. C. Gossard, Two-dimensional magnetotransport in the extreme quantum limit. Phys. Rev. Lett. 48, 1559 (1982).
Thouless1982 D. J. Thouless et al., Quantized Hall Conductance in a Two-Dimensional Periodic Potential. Phys. Rev. Lett. 49, 405 (1982).
Wen1990 X.-G. Wen, Topological orders in rigid states. Int. J. of Mod. Phys. B 4, 239-271 (1990).
Landau1999 L. D. Landau, E. M. Lifshitz, and M. Pitaevskii, 1999, Statistical
Physics (Butterworth-Heinemann, New York).

Hasan2010 M. Z. Hasan and C. L. Kane, Topological insulators. Rev. Mod. Phys. 82, 3045 (2010).
Qi2011 X. L. Qi and S. C. Zhang, Topological insulators and superconductors. Rev. Mod. Phys. 83, 1057 (2011).
Konig2007 M. König et al., Quantum Spin Hall Insulator State in HgTe Quantum Wells. Science 318, 766 (2007).
Hsieh2008 D. Hsieh et al., A topological Dirac insulator in a quantum spin Hall phase (experimental realization of a 3D Topological Insulator). Nature 452, 970 (2008).
Xia2009 Y. Xia et al., Observation of a large-gap topological-insulator class with a single Dirac cone on the surface. Nat. Phys. 5, 398 (2009).
Chang2013 C.-Z. Chang et al., Experimental Observation of the Quantum Anomalous Hall Effect in a Magnetic Topological Insulator.  Science 340, 167 (2013).
Xu2015 S.-Y. Xu et al., Discovery of a Weyl Fermion Semimetal and Topological Fermi Arcs. Science 349, 613 (2015).
Lv2015 B. Q. Lv et al., Experimental discovery of weyl semimetal TaAs. Phys. Rev. X. 5, 31013 (2015).

Benevig2017 B. Bradlyn et al., 
Topological quantum chemistry. Nature 547, 298 (2017).

Kitaev2001 A. Kitaev, Unpaired Majorana fermions in quantum wires. Physics-Uspekhi 44, 131 (2001).
Reed-Green2000 N. Read and D. Green, Paired states of fermions in two dimensions with breaking of parity and time-
reversal symmetries and the fractional quantum Hall effect. Phys. Rev. B 61, 10267 (2000).

Alicea2012 J. Alicea, New directions in the pursuit of Majorana fermions in solid state systems. Rep. Prog. Phys. 75, 76501 (2012).

Elliott2014 S. R. Elliott, M. Franz, Colloquium: Majorana Fermions in nuclear, particle and solid-state physics. Rev. Mod. Phys. 87, 137 (2014).

Sato2017 M. Sato and Y. Ando, Topological superconductors: a
review. Rep. Prog. Phys. 80, 076501 (2017).

Chen2012 X. Chen et al., Symmetry-Protected Topological Orders in Interacting Bosonic Systems. Science 338, 1604 (2012).
Chen2013 X. Chen et al., Symmetry protected topological orders and the group cohomology of their symmetry group. Phys. Rev. B 87, 155114 (2013).

quench1  L. D'Alessio and M. Rigol, Dynamical preparation of Floquet Chern insulators. Nat. Commun. 6, 8336 (2015).

quench2  M. D. Caio, N. R. Cooper, and M. J. Bhaseen, Quantum Quenches in Chern Insulators. Phys. Rev. Lett. 115, 236403 (2015).

quench3  Y. Hu, P. Zoller, and J. C. Budich, Dynamical Buildup of a Quantized Hall Response from Nontopological States. Phys. Rev. Lett. 117, 126803 (2016).

quench4  F. N. Ünal, E. J. Mueller, and M. Ö. Oktel, Nonequilibrium fractional Hall response after a topological quench. Phys. Rev. A 94, 053604 (2016).

quench5  J. H. Wilson, J. C.W. Song, and G. Refael, Remnant Geometric Hall Response in a Quantum Quench. Phys. Rev. Lett. 117, 235302 (2016).


Su1980 W. P. Su et al., Soliton excitations in polyacetylene. Phys. Rev. B 22, 2099 (1980).

Atala2012 M. Atala et al., Direct Measurement of the Zak phase in Topological Bloch Bands. Nat. Phys. 9, 795 (2012).

Liu2013 X.-J. Liu, Z.-X. Liu, and M. Cheng, Manipulating Topological Edge Spins in a One-Dimensional Optical Lattice. Phys. Rev. Lett. 110, 076401 (2013).

Song2017 B. Song, L. Zhang, C. He, T. F. J. Poon, E. Hajiyev, S. Zhang, X.-J. Liu, and G.-B. Jo,
Observation of symmetry-protected topological band with ultracold fermions. Science Advances 4,
eaao4748 (2018).

Jotzu2014 G. Jotzu et al., Experimental realization of the topological Haldane model with ultracold fermions. Nature 515, 237-240 (2014).
Aidelsburger2015 M. Aidelsburger et al., Measuring the Chern number of Hofstadter bands with ultracold bosonic atoms. Nat. Phys. 11, 162-166 (2015).

Wu2016 Z. Wu, L. Zhang, W. Sun, X.-T. Xu, B.-Z. Wang, S.-C. Ji, Y. Deng, S. Chen, X.-J. Liu, and J.-W. Pan, Realization of two-dimensional spin-orbit coupling for Bose-Einstein condensates.  Science 354, 82 (2016).

Sun2017 W. Sun et al., Long-lived 2D Spin-Orbit coupled Topological Bose Gas. arXiv:1710.00717.

Tarnowski2017 M. Tarnowski et al., Characterizing topology by dynamics: Chern number from linking number. arXiv: 1709.01046.

Wang2017 C. Wang et al., Scheme to Measure the Topological Number of a Chern Insulator from Quench Dynamics.  Phys. Rev. Lett. 118, 185701 (2017).

Sunwei2018 W. Sun et al., Uncovering topology by quenching dynamics. arXiv:1804.08226v2.

AZ1997 A. Altland and Martin R. Zirnbauer, Nonstandard symmetry classes in mesoscopic normal-superconducting hybrid structures. Phys. Rev. B 55, 1142 (1997).
Schnyder2008 A. P. Schnyder et al., Classification of topological insulators and superconductors in three spatial dimensions. Phys. Rev. B 78, 195125 (2008).
Kitaev2009 A. Kitaev, Periodic table for topological insulators and superconductors. AIP Conference Proceedings 1134, 22 (2009).
Chiu2016 C.-K. Chiu et al., Classification of topological quantum matter with symmetries. Rev. Mod. Phys. 88, 035005 (2016).
Morimoto2013 T. Morimoto et al., Topological classification with additional symmetries from Clifford algebras. Phys. Rev. B 88, 125129 (2013).
Chiu2013 C.-K. Chiu et al., Classification of topological insulators and superconductors in the presence of reflection symmetry. Phys. Rev. B 88, 075142 (2013).

SI See Supplementary Material for more details.

Haldane1988 F. D. M. Haldane, Model for a Quantum Hall Effect without Landau Levels: Condensed-Matter Realization of the “Parity Anomaly”. Phys. Rev. Lett. 61, 2015 (1988).
Volovik2003 G. Volovik, The universe in a helium droplet (Clarendon Press, 2003).

Zhang2001 S. C. Zhang and J. Hu, A four-dimensional generalization of the quantum Hall effect. Science, 294, 823-828 (2001).



Liu2014 X.-J. Liu, K. T. Law, and T. K. Ng, Realization of 2D
Spin-Orbit Interaction and Exotic Topological Orders in Cold Atoms. Phys. Rev. Lett. 112, 086401 (2014).

Wang2018 B.-Z. Wang, Y.-H. Lu, W. Sun, S. Chen, Y. Deng, X.-J. Liu, Dirac-, Rashba-, and Weyl-type spin-orbit couplings: Toward experimental realization in ultracold atoms. Phys. Rev. A bf 97, 011605(R) (2018).

Lindblad G. Lindblad, On the generators of quantum dynamical semigroups. Commun. Math. Phys.
48, 119-130 (1976).




Shiozaki2014 K. Shiozaki et al., Topology of crystalline insulators and superconductors. Phys. Rev. B 90, 165114 (2014).

Slager2013 R. -J. Slager 
et al., The space group classification of topological
band-insulators. Nat. Phys. 9, 98 (2013).

Chan2017-2 C. Chan, L. Zhang, T.-F. J. Poon, Y.-P. He, Y.-Q. Wang, and X.-J. Liu, Generic Theory for Majorana Zero Modes in 2D Superconductors. Phys. Rev. Lett. 119, 047001 (2017).









§ SUPPLEMENTARY MATERIAL:
NON-EQUILIBRIUM CLASSIFICATION OF TOPOLOGICAL QUANTUM PHASES


  
 Sequation
 Sfigure
 
 In this Supplementary Material, we shall show the elementary representation matrices of Clifford algebra, and provide the detailed proof of bulk-surface duality in the simple cases, as well as more details of numerical results.



§ THE ELEMENTARY REPRESENTATION MATRICES OF CLIFFORD ALGEBRA


The elementary representation matrices of Clifford algebra {γ_i,γ_j}=2δ_ij1 are of dimension n_d=2^d/2 (or 2^(d+1)/2) if the dimensionality of the system d is even (or odd). Here we give an example to construct these Clifford matrices and discuss their trace property, as shown in Tab. <ref>.






§ PROOF OF BULK-SURFACE DUALITY IN THE SIMPLE CASES

To better understand the proof for the generic bulk-surface duality for any dimensions in the Appendix, in this section, we provide the details of the proof for the simple cases, namely, the 1D, 2D and 3D topological phases.

Below, the proof is directly based on the deformed Hamiltonian

    ℋ( k)=h⃗'( k)·γ⃗=sinΘ_ kγ_0+cosΘ_ k∑^d_i=1ĥ_so,i( k)γ_i,

which is topologically equivalent to the basic Hamiltonian ℋ( k)=h⃗( k)·γ⃗ considered in the main text. Without loss of the generality, h_0 defines the band-inversion surfaces (BISs) and h_so≡(h_1,…,h_d) is the spin-orbit (SO) field. Here, as discussed in the Appendix, Θ_𝐤∈[-π/2,π/2] is a smooth and
monotonic function of h_0(𝐤), which approaches the step function and maps
h_0(𝐤) to ±π/2 if h_0(𝐤)≷0
or to 0 if h_0(𝐤)=0, see Fig. 5d in the main text.
For example, it can be chosen as the error function πerf(ξ h_0(𝐤))/2
with ξ→∞.
Therefore, this deformed Hamiltonian extracts all important informations on the BISs, including the normalized SO field ĥ_so=h_so/|h_so|.
Note that the deformed Hamiltonian has been normalized. Accordingly, the winding number or Chern number can be defined using the normalized Hamiltonian <cit.>.



 §.§ 1D topological phases characterized by the winding number


We start with the 1D topological phases, which are the simplest cases, with the basic Hamiltonian
for the topology being two-band models. Accordingly, the Clifford matrices γ_i
are generated by the Pauli matrices, γ_i=0,1∈{σ_x,σ_y,σ_z}.
The Hamiltonian ℋ(k) has chiral symmetry γ=γ_0γ_1 and is characterized by the 1D winding number:

    ν_1   ≡   1/4π∫_ BZTr[γℋℋ]
       =   1/4π∫_BZ k (h'_0∂_kh'_1-h'_1∂_kh'_0)Tr[γγ_0γ_1],

here we have used the fact ∑_i=0^dh_i^'2=1 for
the normalized Hamiltonian (<ref>). With the deformed Hamiltonian and
the trace property of Clifford matrices Tr[γγ_iγ_j]=-2ϵ^ij, the winding number
can be further written as

    ν_1   =   1/2π∫_BZΘ_k ĥ_so
       =   1/2∑_BIS_j[sgn(h_1,R_j)-sgn(h_1,L_j)],

where the unit SO field ĥ_so=sgn(h_1) for 1D, with sgn(x) being the sign function, and the BISs become discrete points (see Fig. 5a in the main text), with the left hand and right hand points of the j-th BIS being denoted as L_j and R_j, respectively. Thus the right hand side of the above formula denotes a zero-dimensional invariant (0-th Chern number). Note that the region 𝒱_BIS enclosed by each
BIS refers to that with h_0(k)<0, see Fig. 5c in the main text. The integral of Θ_k
must become +π at the right-hand (R_j) and -π at the left-hand (L_j) point of
the j-th BIS, which results in the prefactors of sgn(h_1,R_j) and sgn(h_1,L_j), respectively.
For the case with only a single BIS, the equation (<ref>) gives the result with a simple picture that the system is topologically nontrivial
when the SO field has opposite signs at the R and L points of BIS, otherwise it is trivial.

As an example, we consider the 1D model of the AIII class in the main text, which has the Hamiltonian ℋ(k)=h⃗(k)·σ⃗ with h_x=t_sosin k≡ h_so, h_y=0 and h_z=m_z-t_0cos k≡ h_0. Here t_0>0 and
t_so>0 denote the nearest-neighbor spin-conserved and spin-flipped hopping coefficients, m_z is the effective magnetization. When |m_z|>t_0, this model is trivial as there is no BISs. On the other hand, the BIS is determined by k_BIS=±arccos m_z/t_0 when |m_z|<t_0. It is easy to check that h_so>0 for the right hand point of the BIS and h_so<0 for the left hand point, which indicates the winding number ν_1=1.



 §.§ 2D topological phases characterized by the first Chern number


For the 2D topological phases, a two-band model can be formulated with three Pauli matrices that γ_i∈{σ_x,σ_y,σ_z}
for i=0,1,2. Accordingly, the topology of the phase is characterized by the first
Chern number

    Ch_1=-/16π∫_BZTr[ℋ(ℋ)^2].

The trace is taken in the spin space. Using the relation ℋ(ℋ)^2=∑_i∑_j≠ kh'_i h'_j∧ h'_kγ_iγ_jγ_k
and the trace property Tr[γ_iγ_jγ_k]=-2ϵ^ijk,
where the Levi-Civita tensor and the wedge are both anti-symmetric, we further obtain the Chern number that

    Ch_1=-1/4π∫_BZ^2𝐤·(h'_0∇ h'_1×∇ h'_2-h'_1∇ h'_0×∇ h'_2+h'_2∇ h'_0×∇ h'_1).

From the deformed Hamiltonian (<ref>) we obtain that

    h'_0∇ h'_1×∇ h'_2    =sinΘ_𝐤cos^2Θ_𝐤∇ĥ_so,1×∇ĥ_so,2
       -sin^2Θ_𝐤cosΘ_𝐤∇Θ_𝐤×(ĥ_so,1∇ĥ_so,2-ĥ_so,2∇ĥ_so,1),
    
    h'_i∇ h'_0×∇ h'_j    =cos^3Θ_𝐤∇Θ_𝐤×ĥ_so,i∇ĥ_so,j.

With the above results the Chern number is recast into

    Ch_1=-1/4π∫_BZ^2𝐤·[sinΘ_𝐤cos^2Θ_𝐤∇ĥ_so,1×∇ĥ_so,2-cosΘ_𝐤∇Θ_𝐤×(ĥ_so,1∇ĥ_so,2-ĥ_so,2∇ĥ_so,1)].

To simplify the above formula, we can choose a special configuration of the Θ_ k function, while the final result is independent of the choice of Θ_ k. For this we let the SO field (∝cosΘ_𝐤) for the deformed Hamiltonian be truncated to the BISs (sinΘ_𝐤=0) by taking Θ_𝐤 to
approach the step function. Then the contribution from the first term vanishes. On the other hand, ∇Θ_𝐤 tends
to a delta function, and the second term will contribute to the Chern
number. Therefore, the whole integral can be reduced to that for the
second term over the region ℛ with non-zero ∇Θ_𝐤 (see Fig. 5c-d in the main text):

    Ch_1    =1/4π∫_ℛ^2𝐤·cosΘ_𝐤∇Θ_𝐤×(ĥ_so,1∇ĥ_so,2-ĥ_so,2∇ĥ_so,1)
        =1/4π∫_-π/2^π/2Θ_k_⊥ cosΘ_k_⊥∑_j∮_BIS_j k_∥ (ẑ×n̂_⊥)·(ĥ_so,1∇ĥ_so,2-ĥ_so,2∇ĥ_so,1)
       =∑_j1/2π∮_BIS_jℓ̂·(ĥ_so,1∇ĥ_so,2-ĥ_so,2∇ĥ_so,1)
       =∑_j1/2πϵ^mn∮_BIS_jℓ̂·ĥ_so,m∇ĥ_so,n,

where we have decomposed the momentum 𝐤
into (k_⊥,k_∥) in the second line, with k_⊥ perpendicular
to the BIS and pointing to the side 𝒱̅_BIS with
h_0(𝐤)>0 (see Fig. 5c in the main text), and m,n=1,2. From the last line of the above formula we know that the final result is independent of the function Θ_ k. Equation
(<ref>) shows that the Chern number in the bulk reduces to the 1D winding of the SO field ĥ_so
on the 1D BISs (band inversion rings). The integral ∑_j1/2πϵ^αβ∮_BIS_jℓ̂·ĥ_so,α∇ĥ_so,β is performed on the 1D BISs, with the positive direction defined so that the closed integral paths denote the boundary of the vector area V_ BIS. The 1D winding number then counts the total vorticity of the SO field at the degenerate points, i.e., the total topological charge of the 2D Dirac points, characterized by ĥ_so=0 within V_ BIS. The above result implies
that, the first Chern number is simply given by the total topological charge
of the SO field in the region enclosed by BISs, namely, in the 𝒱_BIS
with h_0(𝐤)<0.

Here we take the Chern insulator in the main text as an example, whose Hamiltonian can be written as ℋ(𝐤)=h⃗(𝐤)·σ⃗ with h⃗(𝐤)=(t_sosin k_x,t_sosin k_y, m_z-t_0cos k_x-t_0cos k_y). Without loss of the generality, we can decompose the vector field h⃗ into h_0≡ h_z and h_so≡ (h_y,h_x). When m_z>2, there is no BISs and the system is trivial. When 0<m_z<2, a topological charge -1 located at 𝐤=0 in V_ BIS is enclosed by the BIS, and gives the Chern number Ch_1=-1, see Fig. 6a in the main text. Furthermore, another two topological charges +1 at (0,π) and (π,0) are also included when -2<m_z<0, which renders the Chern number being changed to 1+1-1=+1. Finally, when m_z<-2, the region 𝒱_BIS equals the whole BZ, and all the four Dirac points {(0,0),(π,π),(0,π),(π,0)} should be considered and contribute a trivial total topological charge (1+1-1-1=0).

On the other hand, one can also choose that h_0≡ t_sosin k_x, in which case 𝒱_BIS
with h_0(𝐤)<0 corresponds to the region with -π<k_x<0, and then
the SO field is h_so≡(h_z,h_y).
A pair of topological charges with opposite signs, located at 𝐤=(±arccos[(m_z+1)/t_0],π) for -2<m_z<0
and at k=(±arccos[(m_z-1)/t_0],0) for 0<m_z<2, will be created from vacuum,
moves and finally annihilates, as m_z varies from -∞ to
+∞. Only one of the two charges is located in 𝒱_BIS and enclosed by the BIS, giving the Chern number C_1=+1 (for -2<m_z<0) or C_1=-1 (for 0<m_z<2) (see Fig. 6b in the main text). The results are clearly consistent with those given above by choosing h_0=h_z( k).




 §.§ 3D topological phases characterized by the winding number


For 3D topological phases, the minimal Hamiltonian involves four
bands, and the Clifford matrices take the Dirac forms. The Hamiltonian
has chiral symmetry γ=-∏_i=0^3γ_i and is characterized
by the 3D winding number

    ν_3=1/48π^2∫_BZTr[γℋ(ℋ)^3].

Similar to the 2D case, the relation ℋ(ℋ)^3=∑_i∑_j≠ s≠ th'_i h'_j∧ h'_s∧ h'_t
and the trace gives Tr[γγ_iγ_jγ_sγ_t]=-4ϵ^ijst, with which
the winding number takes the following explicit form

    ν_3=-1/2π^2∫_BZ^3𝐤 { h'_0∇ h'_1·(∇ h'_2×∇ h'_3)-∇ h'_0·[h'_1∇ h'_2×∇ h'_3-h'_2∇ h'_1×∇ h'_3+h'_3∇ h'_1×∇ h'_2]}.

From the deformed Hamiltonian (<ref>), the integrand in the winding number can be simplified by

    h'_0∇ h'_1·(∇ h'_2×∇ h'_3)=    sinΘ_𝐤cos^3Θ_𝐤∇ĥ_so,1·(∇ĥ_so,2×∇ĥ_so,3)
        -sin^2Θ_𝐤cos^2Θ_𝐤∇Θ_𝐤·(ĥ_so,1∇ĥ_so,2×∇ĥ_so,3
       -ĥ_so,2∇ĥ_so,1×∇ĥ_so,3+ĥ_so,3∇ĥ_so,1×∇ĥ_so,2),
    
    h'_i∇ h'_0·(∇ h'_j×∇ h'_k)=    cos^4Θ_𝐤∇Θ_𝐤·(ĥ_so,i∇ĥ_so,j×∇ĥ_so,k).







Similar to the derivation done for the 2D phase, we take a special configuration that Θ_𝐤 approaches the step function. The computation of the winding number finally reduces to the integral on the BISs, given by

    ν_3    =1/2π^2∫_ℛ^3𝐤 cos^2Θ_𝐤∇Θ_𝐤·[ĥ_so,1∇ĥ_so,2×∇ĥ_so,3-ĥ_so,2∇ĥ_so,1×∇ĥ_so,3+ĥ_so,3∇ĥ_so,1×∇ĥ_so,2]
        =1/2π^2∫_-π/2^π/2Θ_k_⊥ cos^2Θ_k_⊥∑_j∫_BIS_j^2k_∥·[ĥ_so,1∇ĥ_so,2×∇ĥ_so,3-ĥ_so,2∇ĥ_so,1×∇ĥ_so,3+ĥ_so,3∇ĥ_so,1×∇ĥ_so,2]
        =∑_j1/4π∫_BIS_j^2k_∥·[ĥ_so,1∇ĥ_so,2×∇ĥ_so,3-ĥ_so,2∇ĥ_so,1×∇ĥ_so,3+ĥ_so,3∇ĥ_so,1×∇ĥ_so,2]
        =∑_j1/8πϵ^mnl∫_BIS_j^2k·ĥ_so,m(∇ĥ_so.n×∇ĥ_so,l).

The last line of the above equation denotes the coverage of the SO field ĥ_so( k) over the 2D spherical surface S^2 when k runs over the BISs, which implies that the 3D winding number of the bulk is mapped to the 2D Chern number defined on the 2D BISs.
As is known that the integral 1/8πϵ^mnl∫_BIS_j^2k·ĥ_so,m(∇ĥ_so.n×∇ĥ_so,l)
counts the topological charge of the Weyl points at h_so(𝐤)=0 of the SO field in the vector area V_ BIS enclosed by the BISs.

We also give an example for the 3D topological phases. Consider the
B phase of superfluid ^3He <cit.>, which can be treated as a 3D topological
superconducting phase in class DIII. The Hamiltonian reads ℋ(𝐤)=h_0(𝐤)σ_z⊗1+h_1(𝐤)σ_x⊗σ_z+h_2(𝐤)σ_y⊗1+h_3(𝐤)σ_x⊗σ_x
with h⃗(𝐤)=(𝐤^2/2m-μ,k_x,-k_y,-k_z).
It is easy to see that the system is trivial for μ<0. For the case with μ>0, a single Weyl
point with h_so=(h_1,h_2,h_3)=0 at k=0 and of charge +1 is enclosed by the BIS, which renders the system being topologically nontrivial with winding number ν_3=+1.

With the above detailed proof for the 1D, 2D and 3D topological phases, the generic proof for bulk-surface duality in the Appendix can be better understood.



§ CHARACTERIZING QUANTUM ANOMALOUS HALL INSULATORS




The 2D quantum anomalous Hall (QAH) model generally reads H( k)=h⃗( k)·σ⃗ with h⃗( k)=(m_x+t_ so^xsin k_x, m_y+t_ so^ysin k_y, m_z-t_0cos k_x-t_0cos k_y), which has been realized in Ref. <cit.>.
Here we identify the topological phases by m_z (m_x=m_y=0): (i) 0<m_z<2t_0 with the Chern number C_1=-1; (ii) -2t_0<m_z<0 with C_1=1; (iii) being trivial for |m_z|>2t_0.
The flexibility in the decomposition of the vector h⃗( k) allows us to consider different quench ways to detect the topology.
As discussed in the main text, we can take h_0≡ h_z,  h_ so≡(h_x,h_y), and the quench is performed by varying m_z, which corresponds to quenching the
h_z axis. Alternatively, we can also take h_0≡ h_y and quench the h_y axis by ramping m_y quickly from |m_y|≫ t_0 to 0 while fixing m_z.
These two quench ways leads to distinct observations, but can be both used to distinguish different topological phases.
The only thing to be noted is that the ratios t_ so^x,y/t_0 should be properly adjusted to clarify the observation in the two quench ways.



When quenching the h_z axis, we set t_ so^x=t_ so^y=0.2t_0, and suddenly change m_z from the deep trivial phase |m_z|≫2t_0
to one of the two topological phases.
In Fig. 2a-b of the main text, we have shown the results of quenching to the phase (i). Here we provide the results of the phase (ii) in Fig. <ref>a-b.  As discussed in the former
sections, the vanishing time-averaged spin texture ⟨σ⃗⟩=0 characterizes the BISs, and the dynamical spin-texture field, defined as g(𝐤)=-∂_k_⊥⟨σ⃗⟩ across the BISs, signifies the SO field ĥ_so(𝐤)
and can be used to classify the topology.
Compared to Fig. 2a-b in the main text, the BIS in Fig. <ref>a becomes a ring centered at 𝐤=(π,π) rather than at (0,0). Here the region enclosed by the BIS, i.e. V_ BIS, contains three topological charges of the SO field: one charge -1 located at (0,0) and other two charges +1 at (0,π) and (π,0), which indicates the Chern number C_1=1 and can be reflected by the winding of the spin-texture field g(𝐤) along the BIS.

When quenching the h_y axis, we set t_ so^x=t_0 and t_ so^y=2t_0.
In this quench way, the BISs are the same for different topological regimes, and are identified to two lines k_y=-π and k_y=0.
However the dynamical spin-texture field on the BISs should be different, which classifies the different topology, as shown in Fig. 2c-d for C_1=-1 and in Fig. <ref>c-d for C_1=1. In Fig. 2c-d, the emergent dynamical spin-texture field g(𝐤) winds clockwise along k_y=0 in the σ_x-σ_z plane while the winding along k_y=-π is zero. The situation becomes reverse in Fig. <ref>c-d, where the spin-texture field g(𝐤) winds anti-clockwise along k_y=-π and becomes trivial along k_y=0, which characterizes the Chern number C_1=1.



§ 2D TOPOLOGICALLY TRIVIAL PHASE


It seems that the post-quench phase always has non-trivial topology as long as the BIS is observed. However, that is not true.
As a counterexample, we consider the topologically trivial model H_ tr=h⃗( k)·σ⃗, where the vector field reads h⃗( k)=(t_ sosin k_x,t_ sosin k_y,m_z-t_0cos 2k_x-t_0cos 2k_y).
After the sudden quench from m_z≫ 2t_0 to |m_z|<2t_0, the time-averaged spin textures also exhibit four rings (see Fig. <ref>a), but
the dynamical spin-texture field contribute no winding, charactering the trivial topology (Fig. <ref>b).






§ DETECTING TOPOLOGY VIA DISSIPATIVE DYNAMICS




The dynamical classification can be generalized to the dissipative quantum dynamics at finite temperature. The spin dynamics can be described by the so-called Lindblad mater equation <cit.>. As shown in the main text, after a such quench, the dissipation
drives the occupation of Bloch states to approach the equilibrium result,
and the spin polarization ⟨γ⃗⟩ approaches to be polarized in the reverse direction of the vector h⃗( k) after a long-time evolution.
Therefore, the dissipation can be exploited as
an alternative tool to detect the topology.

In Fig. 4 of the main text, we take the 2D model as an example, and illustrate the detecting method by quenching the h_z axis.
Here we consider another case of quenching the h_y axis (in contrast to the results shown in Fig. 2c-d).
In Fig. <ref>a, we present short-time-averaged spin textures, where two lines L_1,2 that satisfy ⟨σ⃗⟩=0  can be treated as the BISs approximately. After taking long-time averages (Fig. <ref>b), we obtain non-zero spin polarizations ⟨σ_x,z⟩ on the BISs,
which directly reveal the direction of the vector field h⃗( k)= h_ so( k) by h_x,z∝ -⟨σ_x,z⟩. The winding of the SO field on the BISs is readily obtained via the spin polarizations on the two lines L_1,2 (Fig. <ref>c), and the non-zero winding on L_1 or L_2 then characterize the
topological number of the post-quench phase.

99
[S1]Shiozaki2014 K. Shiozaki et al., Topology of crystalline insulators and superconductors. Phys. Rev. B 90, 165114 (2014).

[S2]Chiu2016 C.-K. Chiu et al., Classification of topological quantum matter with symmetries. Rev. Mod. Phys. 88, 035005 (2016).

[S3]Volovik2003 G. Volovik, The universe in a helium droplet (Clarendon Press, 2003).

[S4]Wu2016 Z. Wu, L. Zhang, W. Sun, X.-T. Xu, B.-Z. Wang, S.-C. Ji, Y. Deng, S. Chen, X.-J. Liu, and J.-W. Pan, Realization of two-dimensional spin-orbit coupling for Bose-Einstein condensates.  Science 354, 82 (2016).

[S5]Sun2017 W. Sun et al., Long-lived 2D Spin-Orbit coupled Topological Bose Gas. arXiv:1710.00717.

[S6]Lindblad G. Lindblad, On the generators of quantum dynamical semigroups. Commun. Math. Phys.
48, 119-130 (1976).

[S7]Hu2016  Y. Hu, P. Zoller, and J. C. Budich, Dynamical buildup of a quantized hall response from
nontopological states. Phys. Rev. Lett. 117, 126803 (2016).



