Full-waveform inversion: matching observed to modeled waveforms

What FWI solves

FWI is a LEAST-SQUARES OPTIMIZATION problem. Given:

• Observed seismic data $d_\text{obs}$ (all traces, all times, all offsets)
• Source wavelet estimate
• A starting velocity model $V_0$

FWI finds the velocity model $V$ that MINIMIZES the misfit between observed and simulated data:

where $d_\text{sim}(V)$ is computed by solving the wave equation (full finite-difference or finite-element simulation) for the current velocity model $V$. The sum is over all sources $s$, all receivers $r$, all time samples $t$.

The problem is massive: at field scale, you have thousands of sources, tens of thousands of receivers, tens of thousands of time samples. $d_\text{obs}$ is a multi-terabyte dataset; forward modeling is expensive (each iteration simulates the wavefield for every source); the velocity model has millions of voxels. This is one of the most computationally demanding inversions in any field of science.

Exercise, watch one inversion work

Figure 8.4.1 runs FWI in your browser, with the engine of the FWI lab further down this section, on a model small enough to watch: a smooth background of the kind tomography provides, and in it two small slow bodies, 160 m apart, that the background does not contain. Five shots are fired along the surface. Each iteration solves the wave equation for every shot, measures the misfit between the modelled and the recorded waveforms, sends the residual back from the receivers through the adjoint, cross-correlates the two wavefields to form the gradient, and steps down it by a distance a line search chooses.

1. At iteration 0 the model is the background, so the residual in (d) is exactly the energy the two bodies scattered. The gradient in (c) already marks where they are, an image made by the wave equation itself.
2. Step to iteration 1. The model in (b) is the gradient of iteration 0 scaled into metres per second, and that one step removes about four fifths of the misfit.
3. Run on to iteration 10 at 4 Hz. The misfit in (e) falls to about 2% of its start while the model error stays near 80%: each body comes back about 106 m wide and 62 m/s slow, against its true 33 m and 300 m/s. Fitting the data is necessary; it does not make the model the earth.
4. Switch to 2.5 Hz. The inversion now fits the data to 1% of the starting misfit with one body where there are two. The shortest wavelength this wavelet carries at the bodies, taken at 2.5 times its peak frequency, is about 330 m, and half of it, 167 m, is more than their 160 m spacing. Ten iterations never part them: what could is the top of the wavelet spectrum, which holds only a few percent of its peak amplitude. Run on with these noise-free data, the inversion first draws two bodies at iteration 63, with the misfit at 0.1% of its start and each body only 38 m/s deep.
5. Switch to 6 Hz. Half the shortest wavelength falls to 70 m, and the bodies part at the first step and come back narrower and deeper. That is the resolution of FWI, about half the shortest wavelength in the data (Virieux and Operto, 2009), set by the frequencies the data hold far more than by how long the inversion runs.

Every run here starts from the right background. Figure 8.4.2 and the lab below, Figure 8.4.3, take that away: from a poor starting model, the frequencies inverted first decide whether FWI reaches the truth at all.

Cycle skipping: the FWI failure mode

FWI descends the gradient: each iteration computes the gradient of the misfit with respect to the velocity model and steps downhill. Descent finds the nearest minimum of the misfit, and that is the true model only if the starting model lies in the true model’s basin.

How wide is the basin? Take one arrival that the starting model predicts late by $\Delta t$. The least-squares misfit compares the two traces sample by sample, so it rises as $\Delta t$ grows, until each predicted peak sits on an observed trough; past that point it falls again, because the predicted wiggle begins to line up with the observed wiggle one cycle later. For a single frequency $f$ the misfit is proportional to $1 - \cos 2\pi f \Delta t$, and its maximum, the edge of the basin, lies at half a period, $\Delta t = T/2 = 1/(2f)$. That is the criterion of Virieux and Operto (2009): the starting model must predict each arrival to within half a period, or the inversion settles in a minimum where every wiggle is matched to the wrong cycle. That is cycle skipping.

In numbers: at a peak frequency of 10 Hz half a period is 50 ms, so the starting model’s traveltimes must be right to within 50 ms on the longest and latest paths, not only on average. A reflection at 2 s arrives 50 ms late when the velocity above it is about 2.4% too slow. At 3 Hz the tolerance is 167 ms, more than three times larger. A real wavelet carries a band of frequencies, and its basin closes a little before half its peak period.

Sources of the starting model: (1) traditional NMO-based velocity analysis from stacked seismic (moderate quality, but fine for shallow sections); (2) reflection tomography (detailed smooth model); (3) regional geologic models (long-wavelength constraints); (4) interpolated well velocities; (5) previous cycle of FWI + manual editing. Modern FWI workflows spend 30-50% of the project time on preparing the starting model.

Figure 8.4.2 measures the basin instead of drawing it. Plate (a) shifts one Ricker wavelet against itself: its basin closes at $0.43/f_p$, 0.86 of half its peak period, $1/(2f_p)$. Plate (b) takes the middle shot of the lab below, through the same Marmousi window, and multiplies every velocity under the water by $1 + \gamma$, holding the water as the lab’s own starting model does; each of its 99 points is a full solution of the wave equation at one of the lab’s three bands. At 2.4 Hz the basin runs from −9% to +12%: past its rims the misfit falls as the model moves further from the truth, and descent walks away. At 1.6 Hz the slow rim is at −13%, and at 1.0 Hz the basin covers the whole sweep, −16% to +16%. A starting model 13% slow, as the lab’s is, therefore cycle-skips at 2.4 Hz alone and walks home from 1.0 Hz: the multiscale remedy of Bunks et al. (1995), and the experiment the lab runs.

Multi-scale FWI: the industry standard defense

Cycle skipping is a FREQUENCY-DEPENDENT problem. At low frequencies (long wavelengths), the half-period tolerance is large, starting models that would cycle-skip at 20 Hz are safely within the capture zone at 5 Hz. This insight is what makes MULTI-SCALE FWI possible:

• Low-frequency stage: filter the data to the lowest usable frequencies (typically 3-5 Hz band). Run FWI. The LONG-WAVELENGTH velocity structure is recovered without cycle skipping even from a poor starting model.
• Progressive refinement: expand the frequency band in stages (5-8 Hz, then 8-15 Hz, then 15-25 Hz, etc.). Each stage starts from the converged output of the previous, which is already close to truth at those frequencies.
• High-frequency stage: the final stage runs at the full data bandwidth, adding the fine-scale detail.

This is WHY modern marine 3D acquisition has pushed for LOWER-FREQUENCY CONTENT. Broadband sources (BroadSeis, IsoMetrix, Broadband Plus) specifically target the 3-8 Hz band that FWI needs for robust multi-scale starting. Land acquisition uses vibrators with low-frequency sweeps or dynamite with high-output sources for similar reasons.

Now solve the wave equation yourself

Everything above this point is a description of FWI. What follows is FWI. Figure 8.4.3 solves the 2D acoustic wave equation with a finite-difference scheme, runs the data residual backward through the exact adjoint of that scheme, forms the gradient by correlating the two wavefields, and descends it. Nothing in it is scripted or precomputed: every plate is the state of a solver running in this page, offline, on your device.

The velocity model is the real thing. It is a 7.2 by 1.8 km window cut from the Marmousi2 P-wave model, the benchmark this section has been describing, read from the original SEG-Y at its native 1.25 m sampling and resampled to a 50 m grid. Those are the actual Marmousi velocities, dipping layers, faults and all.

Four experiments. The first fires one shot into the true model, so you can watch the wavefield refract and scatter and the gather build at the receivers. The second fires the same shot into the starting model, runs the difference between the two gathers backward from the receivers, and shows the gradient that difference makes, before the illumination compensation an update applies: at 2.4 Hz, from a start 13 percent slow, most of its weight asks for slower rock, although the start is too slow almost everywhere. The third and fourth are the point of the section: the same deliberately poor starting model, inverted two ways.

You have already walked this misfit landscape in Figure 8.4.2. Now watch a real inversion fall into it. Experiment 3 inverts the 2.4 Hz band alone: its misfit falls by 39 percent, and the model ends 355 m/s RMS from the truth, worse than the 322 m/s it started from, because the start puts the long-path arrivals 0.205 s late at 2.4 Hz, past the 0.18 s edge of the wavelet’s basin though just inside half a period, 0.208 s, and the inversion fits them to the wrong cycles. Experiment 4 inverts 1.0, then 1.6, then 2.4 Hz from the identical start and ends at 156 m/s. Both runs drive the misfit down; only one of them moves the model toward the truth.

The observed data are modelled with the same solver that inverts them. That is the inverse crime (Wirgin 2004), and it flatters any inversion. Two toggles weaken it, though neither removes it: the grid and the physics stay those of the inversion. With noise of 5 percent of each gather's peak, the misfit barely falls, by under 1 percent in every band, because no model can fit noise; yet the multi-scale model still ends at 158 m/s, so a misfit that hardly moves is no proof of failure either. With a wavelet whose peak frequency is 45 percent too high, every band's misfit falls by about 1 percent and the multi-scale model ends at 324 m/s, slightly worse than its start: the inversion cannot tell a wrong wavelet from a wrong model.

When FWI is worth running

• Complex structural settings: salt basins, sub-salt targets, sub-thrust plays. Conventional tomography cannot resolve the complex velocity structure. FWI delivers detailed velocity models that feed into accurate depth migration.
• Near-surface problems: unconsolidated weathering zones, karst, permafrost. FWI can image the shallow complexity that degrades deeper imaging.
• High-value targets where imaging matters: billion-barrel subsalt prospects; CO₂ storage pilot projects where precise velocity models enable quantitative monitoring.
• Velocity-model refinement for QI: traditional tomography gives smooth velocity; FWI adds detail that improves pre-stack migration and downstream QI inversion.

When FWI may not be worth it: (1) simple basin geometry where tomography + kirchhoff migration already gives good results; (2) very noisy data where FWI amplifies noise more than signal; (3) budget-constrained projects where the compute cost exceeds the imaging value. For most modern large projects, FWI is now the default, the question is how many iterations and what frequency bandwidth, not whether to run it.

FWI variants and extensions

• Acoustic vs elastic FWI: acoustic assumes only P-waves (simpler, faster, used for velocity-model building). Elastic FWI models P and S waves plus density (more accurate, expensive, used for QI-grade outputs).
• Anisotropic FWI: includes Thomsen ε, δ, γ in the forward model. Essential in basins with VTI shales or HTI fractured reservoirs (Section 8.3).
• Viscoacoustic / viscoelastic FWI: includes attenuation (Q). Important in shallow gas-bearing sediments where attenuation is severe.
• Envelope-based FWI: minimizes the misfit of the wavefield ENVELOPE rather than the waveform itself. Less cycle-skipping-prone; used as a starter for traditional FWI.
• Optimal-transport FWI: uses optimal-transport distance between waveforms instead of L2 norm. Highly robust to cycle skipping but more expensive.
• Time-domain vs frequency-domain: time-domain is more flexible for complex geology; frequency-domain is efficient for narrowband sequential inversion. Modern FWI uses time-domain with frequency selection.