Floquet engineering of electric polarization with two-frequency drive

First, we introduce Floquet theory. We consider a quantum system driven by a time-periodic external field with a Hamiltonian satisfying

$\displaystyle H(t+T)=H(t),$

(2)

where $T=2\pi/\omega$ is the period of external fields. Due to the discrete time-translation symmetry of the system, the solution of the time-dependent Schrödinger equation $i\derivative{t}\ket{\psi(t)}=H(t)\ket{\psi(t)}$ can be taken as an eigenstate of the discrete translation (we set $\hbar=1$). Namely, the general solution of the time-dependent Schrödinger equation can be expanded by the Floquet states $\{\ket{\psi_{\alpha}(t)}\}$ of the form

$\displaystyle\ket{\psi_{\alpha}(t)}=e^{-i\epsilon_{\alpha}t}\ket{u_{\alpha}(t)},\ \ket{u_{\alpha}(t+T)}=\ket{u_{\alpha}(t)},$

(3)

where $\epsilon_{\alpha}$ is called Floquet quasi-energy. Equation (3) shows that Floquet theory is based on an analog of Bloch’s theorem, which is widely known in solid state physics. Namely, discrete spatial-translation symmetry in Bloch’s theorem corresponds to discrete time-translation symmetry in Floquet theory.

We can reinterpret the time-dependent Schrödinger equation in the original Hilbert space $\mathcal{H}$ as a static eigenvalue problem in the extended Hilbert space spanned with the original Hilbert space and an additional Floquet index, by expanding the time-periodic part $\ket{u_{\alpha}(t)}$ of the solution of the Schrödinger equation as a Fourier series. Substituting Floquet states $\ket{\psi_{\alpha}(t)}$ to the time-dependent Schrödinger equation and performing some manipulations, we arrive at

$\displaystyle\sum_{m\in\mathbb{Z}}(H_{n-m}-m\hbar\omega\delta_{m,n})\ket{u_{\alpha}^{n}}=\epsilon_{\alpha}\ket{u_{\alpha}^{n}},$

(4)

where $H_{n}=\frac{1}{T}\int_{0}^{T}\differential te^{in\omega t}H(t)$ and $\ket{u_{\alpha}^{n}}=\frac{1}{T}\int_{0}^{T}\differential te^{in\omega t}\ket{u_{\alpha}(t)}$ are the $n$-th Fourier coefficients of $H(t)$ and $\ket{u_{\alpha}(t)}$, respectively. Equation (4) is an eigenvalue equation in the extended Hilbert space Shirley1965 ; Sambe1973 .

Because $\ket{u_{\alpha^{\prime}}(t)}=\ket{u_{\alpha}(t)}e^{im\omega t}$ also satisfies $\ket{u_{\alpha^{\prime}}(t)}=\ket{u_{\alpha^{\prime}}(t+T)}$, there are multiple representations for the same Floquet state $\ket{\psi_{\alpha}(t)}$ in the extended Hilbert space. This implies that the eigenvalue problem (4) includes many redundant solutions with $\epsilon_{\alpha^{\prime}}=\epsilon_{\alpha}+m\omega$. This redundancy can be removed by block-diagonalizing the matrix in Eq. (4), after which the matrix has the form $(H_{\text{eff}}-m\hbar\omega)\delta_{m,n}$ with the block matrix $H_{\text{eff}}$ being an $m$-independent operator on $\mathcal{H}$.

One approach to perform block-diagonalization and remove such redundancy of a Floquet Hamiltonian is known as van Vleck’s degenerate perturbation theory, which is applicable when the driving frequency of the external fields $\omega$ is much larger than the energy scale of the system Eckardt2015 ; Mikami16 . In such a situation, application of van Vleck’s degenerate perturbation theory leads to the effective Hamiltonian $H_{\textrm{eff}}$ in the form,

$\displaystyle H_{\textrm{eff}}^{\textrm{vV}}=H_{0}+\sum_{m\neq 0}\left(\frac{[H_{-m},H_{m}]}{2m\omega}+\frac{[[H_{-m},H_{0}],H_{m}]}{2m^{2}\omega^{2}}+\sum_{n\neq 0,m}\frac{[[H_{-m},H_{m-n}],H_{n}]}{3mn\omega^{2}}\right)+\mathcal{O}(\omega^{-3})$

(5)

by taking the $-m\hbar\omega\delta_{m,n}$ term as the unperturbed part and treating the $H_{n-m}$ term as a perturbation.

Let us briefly explain Berry phase theory of electric polarization resta-RMP ; vanderbilt-kingsmith . From elementary adiabatic perturbation theory, electric current from the $n$-th band is obtained as

$\displaystyle{\bf\it j}_{n}(\lambda)=\derivative{{\bf\it P}_{n}(\lambda)}{t}=\frac{ie\dot{\lambda}}{(2\pi)^{d}m}\sum_{l\neq n}\int\differential^{d}{\bf\it k}\frac{\bra{u_{n{\bf\it k}}}{\bf\it p}\ket{u_{l{\bf\it k}}}\bra{u_{l{\bf\it k}}}\partial_{\lambda}\ket{u_{n{\bf\it k}}}}{E_{n{\bf\it k}}-E_{l{\bf\it k}}}+\textrm{c.c.},$

(6)

where $\lambda\in[0,1]$ denotes the adiabatic parameter and $\ket{u_{n{\bf\it k}}}$ is the $n$-th Bloch’s state. After integration with respect to $\lambda$, we get Berry phase formula:

$\displaystyle{\bf\it P}(\lambda)=\frac{-ie}{(2\pi)^{d}}\sum_{n:\textrm{occ.}}\int\differential^{d}{\bf\it k}\bra{u_{n{\bf\it k}}}\gradient_{{\bf\it k}}\ket{u_{n{\bf\it k}}},$

(7)

or in the gauge invariant form

$\displaystyle\Delta{\bf\it P}={\bf\it P}(1)-{\bf\it P}(0)=\frac{ie}{(2\pi)^{d}}\sum_{n:\textrm{occ.}}\int\differential^{d}{\bf\it k}\int_{0}^{1}\differential\lambda\innerproduct{\gradient_{{\bf\it k}}u_{n{\bf\it k}}}{\partial_{\lambda}u_{n{\bf\it k}}}+\textrm{c.c.}$

(8)

The integrand function of Eq. (8) is the Berry curvature in the parameter space $({\bf\it k},\lambda)$. Hence, electric polarization is expressed in the form of the Berry phase.

In systems with certain symmetries, polarization often vanishes. If a system has spatial inversion symmetry $\mathcal{I}$, the inversion operation will result in $P\to-P$, so that $P=0$ is required. Similarly, if a system has $n$-fold symmetry $C_{n}$, $P=0$ is also required. Bicircular light (BCL) used in this study is convenient as an external field to break such symmetries. By choosing an appropriate frequency ratio, BCL can break not only $\mathcal{I}$ symmetry but also $C_{n}$ symmetry.

Bicircular light consists of two circular light (CL) waves with different frequencies and opposite chirality. BCL can be expressed in the form of a gauge field $A(t)=A_{x}(t)+iA_{y}(t)$ as

$\displaystyle A(t)=A_{\textrm{L}}e^{in_{1}\omega t}+A_{\textrm{R}}e^{-in_{2}\omega t+i\theta},$

(9)

where $A_{\textrm{L}(\textrm{R})}$ is the amplitude of the left (right) CL and $n_{1},n_{2}$ are the integers representing the frequencies of the two CL waves. As shown in Fig. 1, BCL waves draw various rose curves with $n_{1}$ and $n_{2}$. The number of leaves on the rose curve is determined by the integer ratio of $n_{1}$ to $n_{2}$. In particular, we focus on the case of $(n_{1},n_{2})=(1,n)$ in this paper, where the gauge field is given by

$\displaystyle A(t)=A_{\textrm{L}}e^{i\omega t}+A_{\textrm{R}}e^{-in\omega t+i\theta}.$

(10)

This gauge field has $(n+1)$-fold symmetry, so if a system has $C_{m}$ symmetry, tuning the frequency ratio $n$ so that $n+1$ is coprime to $m$ can break $C_{m}$ symmetry and can induce electric polarization. The parameter $\theta$ is the phase difference between two CL waves and controlling it can rotate the rose pattern drawn by BCL. Rotation of the rose pattern is expected to cause rotation of the electric polarization direction.

Su-Schrieffer-Heeger (SSH) model is a famous model for polyethylene, and consists of 1D tight-binding model with bond alternation su79 . The SSH model has mirror symmetry with respect to the bond center and hence supports no polarization. We show that subjecting the SSH model to two frequency drive with BCL effectively realize Rice-Mele model rice-mele , as illustrated in Fig. 2(a). Rice-Mele model is a representative model for 1D ferroelectrics which breaks inversion symmetry and supports nonzero electric polarization.

The SSH model is described by the two-band Hamiltonian which is given by

$\displaystyle\mathcal{H}=\sum_{j=1}^{N}[(t_{0}+\delta t_{0})c_{j,A}^{\dagger}c_{j,B}+(t_{0}-\delta t_{0})c_{j,A}^{\dagger}c_{j-1,B}]+\textrm{h.c.},$

(11)

where $c_{j,X}~{}(c_{j,X}^{\dagger})$ is the annihilation (creation) operator for site $X$ of the $j$-th unit cell and $t_{0}$ is the hopping parameter. Fourier transformation of creation and annihilation operators yields Bloch Hamiltonian (we set lattice constant to be $a=2$):

$\displaystyle H(k)=2t_{0}\cos k\sigma_{x}-2\delta t_{0}\sin k\sigma_{y},$

(12)

where $\sigma_{i}$ $(i=x,y,z)$ are Pauli matrices. $H(k)$ has inversion symmetry $\mathcal{I}$, i.e., $\sigma_{x}H(k)\sigma_{x}=H(-k)$, so that the electric polarization does not appear¹¹1More precisely, the electric polarization $P$ satisfies $P=-P$ mod $a$ when the inversion symmetry is preserves. This condition is satisfied either by vanishing polarization $P=0$ or by the half of the lattice constant $P=a/2$. In the latter case, the polarization can be nullified by changing the choice of the unit cell.. We can break $\mathcal{I}$ symmetry by applying a $C_{3}$ symmetric electric field, which can be realized by 3-fold BCL:

$\displaystyle A(t)=\real[A_{0}e^{i\omega t}+A_{0}e^{-2i\omega t+i\theta}]$

(13)

Figure 2(b) illustrates the rose patterns drawn by this gauge field $A(t)$ for $\theta=0,\pi/2$. Application of the 3-fold BCL to the system leads to the time-dependent Hamiltonian (we set $e=1$)

$\displaystyle H(k+A(t))=2t_{0}\cos[k+A(t)]\sigma_{x}-2\delta t_{0}\sin[k+A(t)]\sigma_{y}.$

(14)

From Jacobi–Anger expansion $e^{iz\cos\phi}=\sum_{n\in\mathbb{Z}}i^{n}\mathcal{J}_{n}(z)e^{in\phi},$ ($\mathcal{J}_{n}$: Bessel functions of the first kind), the time-dependent Hamiltonian (14) can be written as

$\displaystyle\begin{split}H(t)=2t_{0}\real&\left[\sum_{n,m}i^{n+m}e^{i(n+2m)\omega t}e^{ik}e^{-im\theta}\mathcal{J}_{n}(A_{0})\mathcal{J}_{m}(A_{0})\right]\sigma_{x}\\ -2\delta t_{0}\imaginary&\left[\sum_{n,m}i^{n+m}e^{i(n+2m)\omega t}e^{ik}e^{-im\theta}\mathcal{J}_{n}(A_{0})\mathcal{J}_{m}(A_{0})\right]\sigma_{y}.\end{split}$

(15)

This leads to the Fourier coefficients $H_{m}$ for any BCL amplitude $A_{0}$. Here, assuming a small amplitude and truncating the higher order Bessel functions, we expand Eq. (15) up to the second order of $A_{0}$, which yields

	$$\displaystyle H_{0}=(2-A_{0}^{2})\biggl{[}t_{0}\cos k\sigma_{x}-\delta t_{0}\sin k\sigma_{y}\biggr{]},$$		(16)
	$$\displaystyle H_{\pm 1}=t_{0}\biggl{[}-A_{0}\sin k-\frac{1}{2}A_{0}^{2}\cos ke^{\pm i\theta}\biggr{]}\sigma_{x}-\delta t_{0}\biggl{[}A_{0}\cos k-\frac{1}{2}A_{0}^{2}\sin ke^{\pm i\theta}\biggr{]}\sigma_{y},$$		(17)
	$$\displaystyle H_{\pm 2}=t_{0}\biggl{[}-A_{0}\sin ke^{\pm i\theta}-\frac{1}{4}A_{0}^{2}\cos k\biggr{]}\sigma_{x}-\delta t_{0}\biggl{[}A_{0}\cos ke^{\pm i\theta}-\frac{1}{4}A_{0}^{2}\sin k\biggr{]}\sigma_{y},$$		(18)
	$$\displaystyle H_{\pm 3}=-\frac{1}{2}t_{0}A_{0}^{2}\cos ke^{\pm i\theta}\sigma_{x}+\frac{1}{2}\delta t_{0}t_{0}A_{0}^{2}\sin ke^{\pm i\theta}\sigma_{y},$$		(19)
	$$\displaystyle H_{\pm 4}=-\frac{1}{4}t_{0}A_{0}^{2}\cos ke^{\pm 2i\theta}\sigma_{x}+\frac{1}{4}\delta t_{0}t_{0}A_{0}^{2}\sin ke^{\pm 2i\theta}\sigma_{y}.$$		(20)

By computing commutators $[H_{-m},H_{m}]$, the effective Hamiltonian [Eq. (5)] up to the first order of $1/\omega$ is obtained as

$\displaystyle H_{\textrm{eff}}=H_{0}+\frac{3}{2\omega}t_{0}\delta t_{0}A_{0}^{3}\sin\theta\sigma_{z}+\mathcal{O}(\omega^{-2},A_{0}^{4}).$

(21)

This shows that a mass term ($\propto\sigma_{z}$) appears and Rice-Mele model is realized effectively by the two frequency drive.

Since a mass gap opens in the effective Hamiltonian, we can define electric polarization of the lower band by using the Berry phase formula as $P(\theta)=-i\int_{\textrm{BZ}}\frac{\differential k}{2\pi}\bra{u_{-}(k,\theta)}\partial_{k}\ket{u_{-}(k,\theta)}$, where $\ket{u_{-}(k,\theta)}$ is the wave function of the lower band. This gives the electric polarization of the driven system if the lower band of the effective Hamiltonian is fully occupied in the steady state under the driving²²2In the nonequilibrium steady state under the driving, the electron distribution generally depends on the relaxation process. For example, with small dissipation, the energy of the system increases due to the energy gain from the drive and its temperature may approach infinity DAlessio2014 ; Lazarides2014 . In order that the lower energy states are (fully) occupied in the steady state, we need efficient relaxation processes such as phonon scattering.. In a two-level system with a Hamiltonian $H={\bf\it R}\dotproduct{\bf\it\sigma}$, Berry connection of the lower band is written as

$\displaystyle a_{-}(k)=\frac{1}{2}(1+\cos\xi)\partialderivative{\eta}{k},$

(22)

where ${\bf\it R}=R(\sin\xi\cos\eta,\sin\xi\sin\eta,\cos\xi)$. In this model, ${\bf\it R}=((2-A_{0}^{2})t_{0}\cos k,-(2-A_{0}^{2})\delta t_{0}\sin k,\frac{3}{2\omega}t_{0}\delta t_{0}A_{0}^{3}\sin\theta)$, which leads to

$\displaystyle P(\theta)=-\frac{3}{8\pi}\frac{1}{\omega}t_{0}\delta t_{0}A_{0}^{3}\sin\theta\int_{-\pi/2}^{\pi/2}\differential k\frac{1}{|{\bf\it R}|}\frac{\sec^{2}k}{1+(\delta t_{0}/t_{0})^{2}\tan^{2}k}.$

(23)

Figure 2(c) plots the electric polarization $P$ as a function of $\theta$, and as Eq. (23) indicates, if $A_{0}\ll 1$, then $P$ is proportional to $-\sin\theta$. Since $P$ is expressed as an odd function of $\theta$, we can dynamically control the sign of the electric polarization by tuning the sign of $\theta$.

Next, we consider a 2D toy model defined on a square lattice as illustrated in Fig 3(a). The unit cell consists of 4 sites forming a square, labelled X, Y, Z, and W in Fig. 3(a). The intra-unit-cell hopping amplitude is $t_{0}+\delta t_{0}$ and the inter-unit-cell hopping amplitude is $t_{0}-\delta t_{0}$. This model can be regarded as an extension of 1D SSH model. The Hamiltonian is given by $\mathcal{H}=\sum_{{\bf\it k}}{\bf\it c_{k}}^{\dagger}H({\bf\it k}){\bf\it c_{k}}$ with annihilation operators ${\bf\it c_{k}}=(c_{X{\bf\it k}},c_{Y{\bf\it k}},c_{Z{\bf\it k}},c_{W{\bf\it k}})^{\top}$ and the Bloch Hamiltonian

$\displaystyle\begin{split}H({\bf\it k})\!=\!2\matrixquantity(0&t_{0}\!\cos\!k_{y}\!-\!i\delta t_{0}\!\sin\!k_{y}&0&t_{0}\!\cos\!k_{x}\!+\!i\delta t_{0}\!\sin\!k_{x}\\ t_{0}\!\cos\!k_{y}\!+\!i\delta t_{0}\!\sin\!k_{y}&0&t_{0}\!\cos\!k_{x}\!+\!i\delta t_{0}\!\sin\!k_{x}&0\\ 0&t_{0}\!\cos\!k_{x}\!-\!i\delta t_{0}\!\sin\!k_{x}&0&t_{0}\!\cos\!k_{y}\!+\!i\delta t_{0}\!\sin\!k_{y}\\ t_{0}\!\cos\!k_{x}\!-\!i\delta t_{0}\!\sin\!k_{x}&0&t_{0}\!\cos\!k_{y}\!-\!i\delta t_{0}\!\sin\!k_{y}&0),\end{split}$

(24)

where we set the lattice constant to be $a=2$. This Bloch Hamiltonian $H(k)$ has a 4-fold rotational symmetry $C_{4}$, i.e.,

where $\eta=\arg(k_{x}+ik_{y})$ and $U_{4}$ is the matrix representing the $C_{4}$ symmetry. In a similar manner as in Sec. 3.1, this symmetry can be broken by applying a 3-fold BCL, $A_{x}+iA_{y}=A_{0}(e^{i\omega t}+e^{-2i\omega t+i\theta})$, because of $\mathrm{gcd}(4,3)=1$. The effective Hamiltonian under the driving is obtained by a similar calculation as in the previous section (for details, see Appendix A):

$\displaystyle H_{\text{eff}}=H_{0}+\frac{3}{2\omega}t_{0}\delta t_{0}A_{0}^{3}\biggl{[}\sin\theta(\sigma_{z}\otimes I)-\cos\theta(\sigma_{z}\otimes\sigma_{z})\biggr{]}+\mathcal{O}(\omega^{-2},A_{0}^{4}),$

(25)

where $H_{0}=\frac{2-A_{0}^{2}}{2}H(k)$ and $\otimes$ denotes tensor product defined as

with $2\crossproduct 2$ matrices $A=(a_{ij})$ and $B$. This effective Hamiltonian includes the terms $\propto\sigma_{z}$ and can be regarded as a 2D extension of the 1D Rice Mele model that appeared in the 1D case (see Eq. (21) for comparison). Therefore, application of a 3-fold BCL can induce electric polarization in the square lattice model.

Again assuming that only the lowest band is occupied (with an efficient relaxation), we can obtain the electric polarization from the Berry phase formula, as shown in Fig. 3(b,c). We can control the direction of the polarization $(P_{x},P_{y})$ with the relative phase of two frequency drive $\theta$. As in Fig. 3(b), the fact that $P_{x}$ $(P_{y})$ is an odd (even) function of $\theta$ can be confirmed from the symmetry argument. The effective Hamiltonian Eq. (25) satisfies

$\displaystyle U^{\dagger}H_{\text{eff}}(-k_{x},k_{y},-\theta)U=H_{\text{eff}}(k_{x},k_{y},\theta),$

$\displaystyle U=\sigma_{x}\otimes\sigma_{x},$

(26)

where the unitary transformation $U$ is the mirror operation on $y$-axis. Thus, it is derived that $P_{x}(-\theta)=-P_{x}(\theta)$ and $P_{y}(-\theta)=P_{y}(\theta)$. Eq. (25) indicates that $P_{x}(\theta)$ and $P_{y}(\theta)$ are proportional to $-\sin\theta$ and $-\cos\theta$, respectively, when $A_{0}\ll 1$, as in 1D SSH model. However, when $A_{0}$ becomes larger, the $P_{x}$-$P_{y}$ curve shows anisotropy and looks like a rounded square as shown in Fig. 3(c), where $P=\sqrt{P_{x}^{2}+P_{y}^{2}}$ is peaked at $\theta=\pm\pi/4,\pm 3\pi/4$. Figure 3(d) schematically illustrates the reason why $P_{x}$-$P_{y}$ curve shows a $C_{4}$ symmetric pattern. In Fig. 3(d), the underlying $C_{4}$ symmetry of the square lattice is depicted as squares and the gauge field of BCL with the $C_{3}$ symmetry as equilateral triangles. (We note that the direction of the triangle is different from the rose curve pattern actually drawn by the gauge field.) When a bisection of the triangle coincides with the diagonal of the square lattice, which happens at $\theta=\pm\pi/4,\pm 3\pi/4$, the diagonal behaves as a mirror plane for the overall structure, which constrains the direction of electric polarization in the diagonal direction. As the direction of polarization rotates by $\pi/2$, the $P_{x}$-$P_{y}$ curve shows a rounded square pattern reflecting the $C_{4}$ symmetry of the square lattice.

From a more precise symmetry argument, we can prove that the electric polarization reflects the symmetry of the square lattice as follows. We find that the time-dependent Bloch Hamiltonian $H(k_{x},k_{y},t,\theta)$ satisfies the relationship,

$\displaystyle U_{4}^{\dagger}H(k_{y},-k_{x},t-\pi/2\omega,\theta+\pi/2)U_{4}=H(k_{x},k_{y},t,\theta),$

(27)

where $U_{4}$ defined in Eq. (3.2) is an operator that rotates the system by $\pi/2$. (For details of the derivation and its generalization, see Appendix B.) This equation indicates that the shift of the relative phase $\theta\to\theta+\pi/2$ and the time translation $t\to t-\pi/2\omega$ results in $(-\pi/2)$-rotation of the system characterized by $(k_{x},k_{y})\to(k_{y},-k_{x})$. This immediately leads to the $C_{4}$ symmetry of the electric polarization $P$ which can be expressed as

$\displaystyle P_{x}(\theta+\pi/2)=P_{y}(\theta),$

$\displaystyle P_{y}(\theta+\pi/2)=-P_{x}(\theta),$

(28)

as observed in Fig. 3(c). Furthermore, this dynamical symmetry is associated with the unitary operator $U_{4}$ that describes the $C_{4}$ symmetry of the square lattice model in the equilibrium [Eq. (3.2)]. This implies that the symmetry of the electric polarization with respect the phase $\theta$ reflects the symmetry of the underlying lattice of the system, which generally holds as we see for the case of $C_{3}$ symmetric lattice model in the next subsection.

$C_{5}$ symmetric BCL $A_{x}+iA_{y}=A_{0}(e^{i\omega t}+e^{-4i\omega t+i\theta})$ can also break $C_{4}$ symmetry of the 2D model, because of $\mathrm{gcd}(4,5)=1$. By using up to the third order of $A_{0}$ in the Fourier coefficients, we can obtain another effective Hamiltonian including mass terms,

$\displaystyle H_{\text{eff}}=H_{0}+\frac{5}{144\omega}t_{0}\delta t_{0}A_{0}^{5}\biggl{[}\sin\theta(\sigma_{z}\otimes I)+\cos\theta(\sigma_{z}\otimes\sigma_{z})\biggr{]}+\mathcal{O}(\omega^{-2},A_{0}^{6}).$

(29)

This is almost the same as the effective Hamiltonian in the case of $C_{3}$ symmetric BCL driving [Eq. (25)]. Figure 3(e,f) shows the results of the electric polarization when $C_{5}$ BCL is driven on the 2D square lattice. Because of the need for higher order of Fourier coefficients $H_{m}$ and the BCL amplitude $A_{0}$ to produce the mass terms, the mass terms are smaller than in the case of $C_{3}$ symmetric BCL for the same $A_{0}$. Thus the electric polarization is also smaller and anisotropy of the electric polarization almost disappears. Note that the rotation direction of the electric polarization is opposite to that of the $C_{3}$ symmetric BCL driven case due to the difference in the shape of the rose patterns drawn by BCLs. Appendix B gives a more description of the relationship between the shape of the rose curve and the rotation direction of the electric polarization.

Finally, we apply our method to Floquet-engineer electric polarization to the honeycomb lattice. Figure  4(a) shows the structure of the honeycomb lattice. Primitive lattice vectors ${\bf\it a}_{1}$ and ${\bf\it a}_{2}$ are defined as

$\displaystyle{\bf\it a}_{1}=\left(\frac{a}{2},\frac{\sqrt{3}a}{2}\right),\ {\bf\it a}_{2}=\left(\frac{a}{2},\frac{-\sqrt{3}a}{2}\right),$

(30)

where $a$ is the lattice constant. Each unit cell, which is grayed out in Fig. 4(a), contains two sites, labeled A and B. We further introduce a potential difference $2U$ between A site and B site. The tight binding Hamiltonian of our honeycomb lattice model is written as where $\gamma$ is a hopping parameter and $f({\bf\it k})$ represents the hopping from site A to the three nearest-neighbor B sites and can be written as follows using the relative distance to the nearest-neighbor sites ${\bf\it\delta}_{1}=(0,a/\sqrt{3}),~{}{\bf\it\delta}_{2}=(a/2,-a/2\sqrt{3}),~{}{\bf\it\delta}_{3}=(-a/2,-a/2\sqrt{3})$:

$\displaystyle f({\bf\it k})=\sum_{i=1}^{3}e^{i{\bf\it k}\dotproduct{\bf\it\delta}_{i}}=e^{ik_{y}a/\sqrt{3}}+2e^{-ik_{y}a/2\sqrt{3}}\cos(k_{x}a/2).$

(31)

Since A and B sites are not equivalent due to the on-site potential difference $U$, this model does not have $C_{6}$ symmetry but $C_{3}$ symmetry. On the basis of the same idea as above, it is expected that the application of BCL with $C_{4}$ symmetry to this system can break $C_{3}$ symmetry of the honeycomb lattice and induce the electric polarization (due to $\mathrm{gcd}(3,4)=1$). Thus we consider application of the BCL with $C_{4}$ symmetry, $A(t)=A_{0}(e^{i\omega t}+e^{-3i\omega t+i\theta})$, to the Hamiltonian Eq. (LABEL:eq:_H_honeycomb) with the minimal coupling $k\to k+A(t)$. In this case, we compute the effective Hamiltonian with high frequency expansion by keeping the terms up to the second order with respect to $1/\omega$ in Eq. (5). Figure 4(b) shows the band structure of this honeycomb model with and without $C_{4}$ symmetric BCL drive at $k_{y}=0$. Under the $C_{4}$ symmetric BCL drive, the band gaps are different at K and K^′ points [${\bf\it k}=(4\pi/3a,0),(-4\pi/3a,0)$] because the electronic states at the K and K^′ points couple to the BCL differently depending on their chiralities. This asymmetry in the band structure and the associated wave functions are the origins for the nonzero electric polarization under the driving.

Figure 4(c,d) shows the result of the calculation of the electric polarization obtained from the effective Hamiltonian as described above. For $A_{0}\ll 1$, $P_{x}(\theta)$ and $P_{y}(\theta)$ are almost proportional to $\sin\theta$ and $-\cos\theta$, respectively, which clearly indicates that one can control the direction of the induced polarization by the relative phase $\theta$ of BCL. When $A_{0}$ becomes larger, the $P_{x}$-$P_{y}$ curve shows a $C_{3}$ symmetric pattern as shown in Fig. 4(d), which reflects the underlying $C_{3}$ symmetry of the lattice structure. The magnitude of the polarization $P$ shows maximum at $\theta=\pm\pi/3,\pm\pi$. Since the results above essentially rely on the $C_{3}$ symmetry of the lattice structure, these results suggest that such Floquet engineering of electric polarization is possible in general two-dimensional materials with $C_{3}$ symmetry, which include hexagonal boron nitride (BN), transition metal dichalcogenides such as MoS₂, and bilayer graphene with layer asymmetry.