Automatic velocity picking by simulated annealing

Simulated annealing algorithm

The simulated annealing (SA) algorithm is the computational analog of slowly cooling a metal so that it adopts a low-energy, crystalline state. It is a provably convergent optimizer. Geman and Geman (1984) provided a proof that simulated annealing, if the annealing is sufficiently slow, converges to the global optimum. In geophysics, SA has been employed to solve the problems of statics (Rothman, 1985), waveform inversion (Sen and Stoffa, 1991) and ray tracing (Bona et al., 2009). Here we propose the application of the simulated annealing algorithm to automatic velocity picking.

We initialize the system at a high temperature. At this stage the particles are then free to move around by a small pertubation; as the temperature is lowered, however, they are increasingly confined due to the high energy cost of movement. At each temperature T, SA perturbs the system randomly until it reaches equilibrium. The new state of the perturbed system is accepted according to the metropolis algorithm:

$$P =  min(1, e^{-\frac{(E(\rho') - E(\rho))}{T}}),$$ (1)

where $\rho$ is the current velocity model, $\rho'$ is the perturbed velocity model, and E is the objective function presented later in the chapter. As shown in the equation, a higher T introduces a higher probability that an uphill perturbation will be accepted, which means that SA can draw samples from the whole population. As T decreases, only perturbations leading to smaller increases in E are accepted, so that only limited exploration is possible as the system settles to the global minimum. The pseudo code of the algorithm is given in Table 1.

Table 1: Algorithm - Simulated Annealing

Inputs:

$\{T_{k}\}_{k=1}^{K}$  Sequence of temperature values

$\{L_{k}\}_{k=1}^{K}$  Sequence of epoch durations

$\rho$         Initial velocity model

SA loop:

   1:    do itemp = 1, K

   2:        do iepoch = 1, $L_{K}$

   3:            $\rho' : = perturb(\rho)$

   4:            $\delta$E: = $E(\rho') - E(\rho)$

   5:            rand: = rand(0,1)

   6:            if rand $<$ min(1, exp(-$\delta$E/$T_k$))

   7:                $\rho$ = $\rho'$

   8:            end if

   9:        end do

  10:     end do

The simulated annealing algorithm is popular and has been well-developed for single-objective optimization. Traditionally, multiobjective optimization can be converted into a single-objective optimization by different fix-up approaches such as the weighted-sum or $\epsilon$-constraint method. In our case, we use a composite-objective function:

$$E(\rho) = \omega_{semb}  \frac{1}{Semb(\rho) + \epsilon} +  \omega_{smooth}  \nabla \rho,$$ (2)

where $Semb(\rho)^{-1}$ is the inverse of semblance and is a measurement of focusing; $\epsilon$ is a arbitrary small number to avoid diverging when $Semb(\rho) = 0$; and $\nabla \rho$ is the residual after passing the velocity model through the Laplacian operator and is a measurement of smoothness. $\omega_{semb}$ and $\omega_{smooth}$ are the weights for these two measurements, respectively, making $E(\rho)$ is the linear combination of these two objective functions.

There are different ways to compose this single objective function. By using the inverse of semblance and residual of Laplacian we cast the optimization as a minimization problem. This combined objective function is then used as the energy or cost to be minimized in an SA optimizer. For the case we discuss here, the fixed weights are chosen according to the relative importance of the objective functions.

Automatic velocity picking by simulated annealing