## 3 Results
### Domain Model  
As is already mentioned, the treatment based on the extended Hubbard model has a disadvantage of a restricted system size. To understand the dynamical behaviors depending on \(\Delta_{eff}\), the use of the domain model [12; 13; 14] is considered to be effective. The domain model is defined as  
\[H_{eff}=-t_{0}\sum_{l_{D},p}\{(|l_{D},p+1\rangle+|l_{D}-1,p+1\rangle)\langle l _{D},p|+h.c.\}+\sum_{l_{D},p}E(p)|l_{D},p\rangle\langle l_{D},p|\, \tag{7}\]  
where \(|l_{D},p\rangle\) is the domain state of which the spacial size and the position of the leftmost site are \(p\) and \(l_{D}\), respectively. The domain energy, \(E(p)\), is defined as \(E(p)=V+p(2V+\Delta_{eff}-U)\), which consists of the domain internal energy and the domain-wall energy. While this model considers only the domain degree of freedom, it is known to reproduce the properties of the extended Hubbard model close to the phase boundary when the problem is limited to the case in which only one domain exists, that is, the limit of weak excitation. The advantage of this model is, of course, the nature of a one-body problem, which makes calculations in large systems realistic. Another remark related to the treatment of this model is the use of a first-order perturbational theory. Since we assume weak excitation as we already mentioned, this treatment is considered to be appropriate. In actual numerical calculations, we slightly modify the form of \(A(t)\) by adding a factor of \(\exp(-\gamma t)\) to avoid any accidental divergence in the energy denominator. We set \(\gamma\) at 0.002 eV, and it is confirmed that the final result is almost the same as the result with \(\gamma\) = 0 when the latter is safely obtained. Furthermore, we restrict our calculation in the subspace of zero total momentum because the momentum of visible or infrared light is negligible. In this case, the effective model is further reduced to that with only \(p\), which is expressed as  
\[H_{eff}=-2t_{0}\sum_{p}((|p+1\rangle\langle p|+h.c.)+\sum_{p}E(p)|p\rangle \langle p|. \tag{8}\]  
Figure 5: (**a**–**c**) Same as Figure 3, but for \(\Delta_{eff}\) = 0.70 eV and \(\hbar\omega_{c}\) = 1.2 eV. We emphasize that the transfer energy is doubled due to the motions of the two domain walls on both sides [14]. Regarding the effect of photoexcitation, it is expressed by the following perturbational term as  
\[H_{e-ph}=E_{0}\,\theta(t)\theta(\tau_{d}-t)\sin(\omega_{c}t)(|p=1)\langle g|+h.c.\rangle. \tag{9}\]  
Since we use the first-order perturbational theory and only discuss the motions in the excited states, the magnitude of electric field \(E_{0}\) is irrelevant in the following. In Figure 6a, we draw the time evolution of the averaged domain size, which is defined as \(\sum_{l_{D},p}p|\langle\Psi(t)|l_{D},p\rangle|^{2}\), with \(|\Psi(t)\rangle\) being the time-dependent state vector. \(\hbar\omega_{c}\) and \(\tau_{d}\) are 1.7 eV and 18 (\(\hbar/\)eV), respectively. As is seen immediately, the oscillatory behavior is observed and its period gradually increases toward the critical point at \(\Delta_{eff}\) = 0.27 eV. From a different viewpoint, these are interpreted as the motions of a wave packet along the relative coordinate of the two domain walls, as schematically depicted in Figure 6b. To understand the nature of these oscillations, we perform Fourier analyses, particularly, for \(\Delta_{eff}\) = 0.50 eV. In Figure 7a, we show the oscillatory components from which the background components are already subtracted. The red curve is obtained from the corresponding curve in Figure 6a, while the blue curve from the result with \(\hbar\omega_{c}\) being as high as 10.2 eV. The comparison of these two cases in Figure 7a does not give any meaningful information. However, their Fourier transforms in Figure 7b show that the numbers of principal peaks are two and one for the former and the latter, respectively.