In order to relieve extensive data storage and burdensome disk I/O and thus reach a reasonable compromise between the computer memory requirement and computational complexity, we propose an efficient wavefield reconstruction strategy named CATRC, which combines the efficiency of reverse propagation and the stability of checkpointing. Therefore source wavefields used in imaging condition in equation 4 can be well-reconstructed during time-reversal simulation. Here we denote the reconstructed wavefields as $q(\mathbf{x},t)$, which is the solution of

...ig \} \bigcup \big \{T, T-\Delta t\big \}, &
\end{array}\right.\end{displaymath} (5)

where $\partial \Omega$ is the boundary of space $\Omega$, $h(\mathbf{x},t)$ are forward wavefields distributed on $\partial \Omega$, which are reversed in time and enforced as the Dirichlet boundary condition for source wavefield reconstruction. $C_i$ denotes a checkpoint, and $N_c$ is the number of checkpoints. According to equation 5, the implementation of CATRC can be briefly summarized as two processes. Firstly, we compute forward wavefield $p_s(\mathbf{x},t)$ by solving the compensated viscoacoustic wave equation in 2, and we save the forward wavefield at the outermost layer boundary of the simulation domain at every time step. At the same time, we also save the complete forward wavefield $p_s(\mathbf{x},t)$ at the predefined checkpoints ( $t \in C_i, i=0, \cdots , N_c-1$) and the last two time steps. The checkpoints can be equally distributed or logarithmically distributed (Symes, 2007; Griewank and Walther, 2000). Next, we compute the backward wavefield $q(\mathbf{x},t)$ in reverse time (from $t=T$ to $t=0$) by solving the reconstructed viscoacoustic wave equation (equation 5), and replacing the calculated backward wavefield $q(\mathbf{x},t)$ with the recorded forward wavefield $p_s(\mathbf{x},t)$ at checkpoints ( $t \in C_i, i=N_c-1, \cdots , 0$).

It is remarkable that reconstruction by equation 5 is a mathematically stable process, given that the source wavefield is compensated while the reconstructed wavefield from boundary is attenuated. However, this stable reconstruction still suffers from insufficient accuracy due to the fact that we utilize PSM to solve equation 5 with only the recorded forward wavefield at the outermost layer boundary of simulation domain at every time step. This mismatch of simulation accuracy inevitably degrades the performance of the wavefield reconstruction. Fortunately, the time-reversal checkpointing scheme acts as a time-domain regularization that eliminates accumulating errors by replacing the reconstructed wavefield with the stored wavefield at checkpoints (Wang et al., 2017b).