--- title: Optimized time-lapse acquisition design via spectral gap ratio minimization author: | Yijun Zhang^1^, Ziyi Yin^1^, Oscar Lopez^2^, Ali Siahkoohi^3^, \ Mathias Louboutin^1^, Rajiv Kumar^4^, Felix J. Herrmann^1^ \ ^1^ Georgia Institute of Technology, ^2^ Florida Atlantic University, ^3^ Rice University, ^4^ Schlumberger bibliography: - abstract.bib --- ## Abstract: Modern-day reservoir management and monitoring of geological carbon storage increasingly call for costly time-lapse seismic data collection. In this letter, we show how techniques from graph theory can be used to optimize acquisition geometries for low-cost sparse 4D seismic. Based on midpoint-offset domain connectivity arguments, the proposed algorithm automatically produces sparse non-replicated time-lapse acquisition geometries that favor wavefield recovery. ## Introduction Time-lapse seismic data acquisition is a costly but crucial endeavor for reservoir management and monitoring of geological carbon storage. While sparse randomized collection of seismic data can lead to major improvements in acquisition productivity [@herrmann2008non;@hennenfent2008simply;@herrmann2010randomized;@mosher2014increasing], systematic approaches to performance prediction, other than computationally expensive simulation-based studies, are mostly lacking. Besides, acquisition optimization approaches, such as minimizing the mutual coherence [@tang2008CITEocs;@mosher2014increasing;@obermeier2017sensing] or minimizing the spectral gap ratio [SGR, @lopez2022graph;@zhang2022simulation], do not handle the unique challenges of time-lapse seismic data acquisition. To meet these challenges, inversion with the joint recovery model [JRM, @oghenekohwo2017low;@wason2017low] will be combined with automatic binary sampling mask generation driven by SGR minimization [@zhang2022simulation]. We opt for the JRM because it inverts baseline and monitor surveys jointly for the common component, which contains information shared between the surveys, and innovations with respect to the common component. Since the fictitious common component is observed by all surveys, its recovery improves when the time-lapse surveys contain complementary information. This is the case when sparse surveys are not replicated [@oghenekohwo2017low;@wason2017low] or when the time-lapse datasets contain independent noise terms [@tian2018joint]. In either case, the JRM leads without insisting on replication of the surveys to high degrees of time-lapse repeatability both in the data [@oghenekohwo2017low;@wason2017low] and image space [@yin2023derisking]. It also yields better interpretability of time-lapse field data [@wei2018improve]. As demonstrated by this letter, including the common component offers additional advantages when optimizing time-lapse acquisition via SGR minimization. To demonstrate this, we first explain the relationship between the SGR and connectivity within graphs associated with binary sampling masks. Next, we describe how this connectivity, which favors wavefield reconstruction, can be improved by minimizing the SGR via optimization. To enhance inversion of time-lapse data with the JRM, a new optimization objective will be introduced that contains SGRs of the common component and of the baseline/monitor surveys. After a brief discussion on minimizing this objective with simulated annealing, the proposed methodology for automatic time-lapse binary mask generation is numerically validated on realistic synthetic 2D data. ## Optimized time-lapse acquisition While the SGR has been used successfully to predict and improve the performance of wavefield reconstruction, it has not yet been used to optimize time-lapse acquisition. After briefly discussing the SGR and JRM, we introduce our methodology to optimize time-lapse data acquisition. ### The spectral gap ratio As shown by @lopez2022graph\, the success of seismic wavefield reconstruction via universal matrix completion [@bhojanapalli2014universal] can be predicted by the ratio of the first two singular values of binary sampling masks, ``{\sigma_2(\mathbf{M})}/{\sigma_1(\mathbf{M})}\in [0,\, 1]`` where ``\mathbf{M}`` is a binary matrix with ``1``'s where data is sampled and with ``0``'s otherwise\. This ratio is known as the spectral gap ratio (SGR) and provides a cheap-to-compute quantitative measure to predict recovery quality. The smaller the SGR, the better the connectivity within graphs spanned by binary sampling masks. Improved connectivity leads to improved wavefield recovery [@lopez2022graph]. While useful, the SGR itself is not constructive because it does not produce sampling masks with small SGRs. With simulated annealing, @zhang2022simulation came up with a practical algorithm to generate acquisition geometries with small SGRs. In this work, we extend this approach by optimizing sparse geometries for time-lapse data acquisition. ### Optimized sampling mask generation Given an initial binary mask, ``\mathbf {M} \in \{0,1\}^{n_s \times n_r}``, with ``n_s`` sources and ``n_r`` receivers, @zhang2022simulation proposed a methodology to minimize the SGR via ```math #SG_opt_1V \underset{\mathbf{M}} {\text{minimize}} \quad \mathcal L(\mathbf{M}) \quad \text{subject to} \quad \mathbf{M}\in \mathcal{C}, ``` with the objective, ``\mathcal L(\mathbf{M}) = {\sigma_2(\mathbf{M})}/{\sigma_1(\mathbf{M})}``, given by the SGR\. To ensure feasibility of the optimized binary masks, the constraint, ``\mathcal{C}=\bigcap\limits_{i=1}^3 \mathcal{C}_i``, is included, which consists of the intersection of the cardinality constraint, ``\mathcal{C}_1=\{\mathbf{M}\mid\#(\mathbf{M})=\lfloor n_s \times \rho\rfloor \times n_r \}``, the binary mask constraint, ``\mathcal{C}_2=\{\mathbf{M}\mid \mathbf{M}\in\{0,1\}^{n_s \times n_r}\}``, and a constraint on the maximum gap size between consecutive samples, ``\mathcal{C}_3=\{\mathbf{M}\mid \operatorname{maxgap}(\mathbf{M})\leq\Delta\}``, where ``\Delta`` is the maximal permitted gap size\. By solving Equation #SG_opt_1V\, @zhang2022simulation produced binary sampling masks that improved wavefield reconstruction compared to masks generated with randomized jittered sampling [@hennenfent2008simply]. Figure #fig:Jittered_vs_Optimized_Mask contrasts jittered with optimized sampling in the midpoint-offset domain, reducing the SGR from ``0.333`` to ``0.196``. The optimized mask increases the sampling at the near offsets where there are more ways to connect to midpoints, which favors wavefield reconstruction [@lopez2022graph]\. #### Figure: {#fig:Jittered_vs_Optimized_Mask} ![](Results/1Vintages/Jittered_MH_mask.png){#J_Mask width=50%} ![](Results/1Vintages/Optimized_MH_mask.png){#P_Mask width=50%} : (a) Jittered versus (b) optimized sampling mask in the midpoint-offset domain. ### Joint recovery model Lowering costs while ensuring time-lapse repeatability are the main challenges in the design of seismic monitoring systems employed to optimize reservoir management and to safeguard geological carbon storage. Both challenges can be met by inverting sparsely sampled baseline and monitor data jointly. For time-lapse acquisition with a single monitor survey, this entails inverting ```math #eq-lin-model \mathbf{b}=\mathcal{A}\left(\mathbf{Z}\right)\quad\text{with}\quad \mathcal{A}\left(\cdot\right) = \begin{bmatrix} \mathcal{A}_1 & \mathcal{A}_1 & 0 \\ \mathcal{A}_{2} & 0 & \mathcal{A}_{2} \end{bmatrix}\left(\cdot\right). ``` In this JRM, the linear operators, ``\mathcal{A}_j,\, j=1, 2``, stand for the combined action of converting monochromatic time-lapse data from the midpoint-offset to the source-receiver domain, followed by trace collection with the acquisition geometries defined by the binary sampling masks, ``\mathbf{M}_j,\, j=1, 2`` with ``j=1`` and ``j=2`` masks for the baseline/monitor surveys. With this model, time-lapse data, ``\mathbf{b}``, which contains the baseline, ``\mathbf{b}_1`` and monitor data, ``\mathbf{b}_2``, are linearly related to ``\mathbf{Z}``, which contains matrices for the unknown densely sampled common component, ``\mathbf{Z}_0``, and innovations with respect to this common component, ``\mathbf{Z}_j,\, j=1, 2``. Compared to recovering the baseline/monitor surveys separately, the JRM produces repeatable results from non-replicated [@oghenekohwo2017low;@wason2017low;@kumar2017highly], non-calibrated [@oghenekohwo2017highly], and noisy [@tian2018joint], time-lapse data. These enhanced results are due to the improved recovery of the fictitious common component. ### Time-lapse optimized mask generation Based on the success of the JRM, we carry the argument of minimizing the SGR a step further by optimizing this quantity for the baseline/monitoring surveys. Because ``\mathbf{Z}_0`` is observed by both surveys, the set of sampling points, ``\left\{\mathbf{M}_0\right\}``, equals the union ``\left\{\mathbf{M}_0\right\}=\left\{\mathbf{M}_1\right\}\cup\left\{\mathbf{M}_2\right\}``. When surveys are replicated, ``\left\{\mathbf{M}_0\right\}=\left\{\mathbf{M}_1\right\}=\left\{\mathbf{M}_2\right\}``. However, ``\mathbf{M}_0`` becomes larger when the baseline and monitor surveys are not replicated explaining why the common component is better resolved when the surveys are not replicated. While Equation #SG_opt_1V leads to improved sampling masks for individual surveys, it does not exploit the fact that the common component is observed by all surveys. For this reason, we propose an optimization procedure with respect to ``\mathbf{M}_1`` and ``\mathbf{M}_2`` with an objective that also includes the SGR for the common component. To avoid generation of poor sampling masks, we follow a mini-max principle where the maximum---i.e., the ``\ell_\infty``-norm---of the SGRs for the common and innovation components is minimized. To compensate for likely smaller SGRs for the common component when the surveys do not overlap (``\#\left\{\mathbf{M}_0\right\}> \#\left\{\mathbf{M}_1\right\}, \#\left\{\mathbf{M}_2\right\}``), we also introduce a scaling. We base this scaling on the property [see Definition ``3.1`` in @bhojanapalli2014universal; @hoory2006expander] that the second singular value of ``d``-regular graphs---i.e., seismic sampling masks with ``d`` non-zero entries per midpoint or offset--- scales with ``\sqrt{d}``. Given this scaling, we propose to minimize the following constrained optimization problem, for ``\, j=1,2``: ```math #SG_opt_2V &\underset{\mathbf{M}_1,\, \mathbf{M}_2} {\text{minimize}} \quad \mathcal{L}({\mathbf{M}}_1, {\mathbf{M}}_2) \quad \text{subject to} \quad \left\{\mathbf{M}_0\right\}=\left\{\mathbf{M}_1\right\}\cup\left\{\mathbf{M}_2\right\}, \mathbf{M}_j\in \mathcal{C}_j,\\ ``` with ``\mathcal L(\mathbf{M}_1, \mathbf{M}_2) = \left\|\left[ \mathcal L(\mathbf{M}_0), \sqrt{\frac{\# (\mathbf{M}_1)}{\# (\mathbf{M}_0)}}\mathcal L(\mathbf{M}_1),\sqrt{\frac{\# (\mathbf{M}_2)}{\# (\mathbf{M}_0)}}\mathcal L(\mathbf{M}_2) \right] \right\|_{\infty}.`` As before, the minimization is subject to constraints, ``\mathcal{C}_j,\, j=1,2``, which can be chosen for each time-lapse survey separately. To produce time-lapse sampling masks, we employ simulated annealing as proposed by @zhang2022simulation but with the following differences: *(i)* randomly perturbed masks are drawn for each survey independently; *(ii)* the compound objectives and constraints of Equation #SG_opt_2V are used; *(iii)* to be relocated sample points are allowed to move more freely than during jitter sampling, which allows us to better explore candidate sampling masks. Figure #fig:intermediate_SG_vs_SNR illustrates how the algorithm progresses when initialized with a replicated jittered subsampled (removing 80% of the sources) acquisition. From Figure #Int_SG\, we observe that the co-located source positions (denoted by the black dots) are gradually replaced by non-coincident source locations for the baseline (blue dots) and monitor surveys (red dots). Even though the objective of Equation #SG_opt_2V decreases non-monotonically (see Figure #SGvsSNR), the reconstruction SNR increases for the baseline and monitor surveys for the selected points. ####Figure: {#fig:intermediate_SG_vs_SNR} ![](Results/2Vintages/intermediate.png){#Int_SG width=58%} ![](Results/2Vintages/SGvsSNR.png){#SGvsSNR width=42%} : Automatic time-lapse sampling mask generation. (a) Starting from a jittered replicated sampling mask, the algorithm produces masks that have smaller SGRs but are no longer replicated. (b) Non-monotonically decaying objective and reconstruction SNR evaluated at points where the objective decreased by more than ``0.003``. ## Numerical validation To confirm the benefits of optimized acquisition, we consider time-lapse data, which differs by a complex gas cloud [@wason2017low;@jones2012building]. Using finite-differences [@witteJUDI2019; @louboutin2022slimgroup; @devito-api; @luporini2020architecture], fully sampled (split spread) 2D baseline and monitor surveys are simulated each consisting of ``300`` sources/receivers sampled at ``12.5 \, \mathrm{m}``. By using a single jittered subsampling mask ``80\%`` of the sources are removed, yielding an average source sampling rate of ``62.5 \, \mathrm{m}`` with ``100\%`` overlap. After running the optimization, the SGRs of the baseline/monitor surveys decreases from ``0.346`` to ``0.268`` and ``0.262``, respectively. The reduction in the overlap ratio (to ``22\%``) leads to improvement in wavefield recovery via matrix completion [@kumar2015efficient; @kumar2017highly], which results in better SNRs for the baseline from ``6.55 \, \mathrm{dB}`` to ``17.03 \, \mathrm{dB}`` and for the monitor from ``6.67 \, \mathrm{dB}`` to ``16.99 \, \mathrm{dB}``. For reasons explained by @oghenekohwo2017low and @wason2017low\, time-lapse difference plots are not included because the benefits of exact replication vanish when acquisition geometries undergo relatively small (``1-2``m) random shifts. ####Figure: {#example_5} ![](Results/2Vintages/jittered.png){#Initial_Jittered_replicated width=50%} ![](Results/2Vintages/Proposed.png){#proposed_nonreplicated width=50%} : Time-lapse wavefield reconstruction in the time-domain. (a) wavefield reconstruction from ``80 \%`` jittered subsampling for the baseline SNR = ``6.55\, \mathrm{dB}``, monitor SNR = ``6.67\, \mathrm{dB}``, and errors between the ground truth and the reconstructed wavefields. (b) the same but with optimized sampling masks, yielding improved recovery baseline/monitor surveys with SNR = ``17.03,\, 16.99\mathrm{dB}``, respectively. While these improvements are encouraging, the proposed optimization is approximate and the produced binary masks will be different for different starting masks. To investigate this effect, ``30`` overlapping jittered masks are generated by removing ``75\%`` of the sources. By reducing the overlap to ``29\% \pm 8\%``, the optimized masks improve the SGRs as can be observed from the violin plots in Figure #SGR\. As before, the reductions in SGRs translate into improved SNRs as can be seen in Figure #SNR\. Compared to box plots, violin plots display the entire distribution including lines for the median (long dashes), first and third quartile (short dashes). We can make the following observations: *(i)* the SGRs for the baseline/monitor surveys decrease significantly; *(ii)* because of the larger number of sampled sources, the SGR for the common component is smaller and more narrowly distributed; *(iii)* the distribution of the SGRs of the baseline/monitor surveys is also narrow compared to the one of the initial jittered binary sampling masks; *(iv)* the SNRs for the recovered baseline/monitor surveys improve significantly. ####Figure: {#fig:different-starting} ![](Results/2Vintages/SG123.png){#SGR width=51%} ![](Results/2Vintages/SNR123.png){#SNR width=49%} : Violin plots for the SGRs (a) and recovery SNRs (b) for ``30`` independent experiments. These experiments show systematic reductions in SGR and significantly improves reconstruction SNRs for optimized surveys. Even though the above results are encouraging and consistent with published reports that claim benefits of the JRM [@oghenekohwo2017low;@wason2017low;@yin2023derisking], further scrutiny is in order. To this end, additional experiments were conducted to better understand robustness of the proposed methodology. Aside from predictable behavior for different starting masks (Figure #fig:different-starting), we also found that optimized SGRs are relatively insensitive to different runs of SA and to random perturbations in the optimized masks. The first observation implies that while SA may produce different masks, the SGRs remain very close, yielding wavefield reconstructions of near equal quality. The second observation indicates that postplot errors by single gridpoint shifts (``12.5``m) in the worst scenario offset the gains made by the optimization. However, on average improvements are mostly preserved although with higher variability. The observed robustness of the presented method is consistent with reported behavior of the JRM. Even though we only considered the on-the-grid case, the argument can be made that improvements will carry over to the off-the-grid situation [@wason2017low;@oghenekohwo2017highly;@lopez2016off]\. However, to turn this claim into a more solid argument, we would have to extend the presented approach to the infinite-dimensional case, which is beyond the scope of this letter. ## Conclusions Acquisition costs form a major impediment to time-lapse seismic. To reduce these costs while ensuring time-lapse repeatability, a novel acquisition optimization scheme was proposed that produces binary sampling masks that favor wavefield reconstruction with the joint recovery model. Optimized sampling masks were generated automatically by minimizing a new objective function consisting of spectral gap ratios for the baseline/monitor surveys and for the common component shared by the surveys. Aside from allowing for wave-simulation free, and therefore computationally feasible, optimized acquisition design, the proposed method also reaffirms the suggestion that deliberate relaxation of survey replication may lead to improved quality of jointly inverted surveys. This claim is solely based on connectivity arguments for the acquisition geometries associated with the baseline/monitor surveys and the common component. Because the spectral gap ratio is extremely cheap to evaluate, it lends itself very well to be extended to multiple monitoring surveys and to 3D. Off-the-grid acquisition geometries are also conducive to being improved by spectral gap ratios, but we will leave that extension to future work. ```math_def \def\argmin{\mathop{\rm arg\,min}} \usepackage{algorithm2e} \usepackage{breqn} \def\minim{\mathop{\hbox{minimize}}} ```