This is an example how to fit the gridded velocities and the gridded interfaces by a model with interfaces, suitable for ray tracing.
We take the gridded velocities and interfaces of the data set 'hess.dat', provided by Amerada Hess, and build the model suitable for ray tracing.
File 'hess.dat' has the form of a file with input
data for the MODEL package. The dimensions of the model are 40200 times
21000 feet. Two interfaces in the model are discretized on the grid of
64 points with spacing 200 feet.
Background velocities are discretized on the grid of
202*106 points with spacing 200 feet*200 feet, velocity inside
a salt body is constant.
Fig. 1: Model 'hess.dat'. |
We can see that the horizontal variation of the velocity field is less pronounced than the vertical variation, and that the vertical variation is greater in the shallow parts of the model than in its deeper parts. The interfaces are quite smooth, but the edges are sharp. According to these observations, we manually choose new parametrization of the model. Let's emphasize, that the parametrization is based on authors opinion, and is expected to be optimized during the inversion. The parametrization is as follows: for the interfaces we take irregular B-spline grid of 28 points with 26 points regularly spaced each 500 feet in the region of the salt body and 2 points at the sides of the model. For the background velocity, we take horizontally regular grid of 17 points with spacing 2500 feet and vertically irregular grid of 31 points with spacing from 200 to 800 feet. The parametrization may be found in file 'hes1-mod.dat', which is used as a starting model for the inversion.
For interfaces we have the data grid of 2*64=128 points with spacing 200 feet and the irregular B-spline grid of 2*28=56 points. For the velocities we have the data grid of 202*106=21412 points with spacing 200 feet. The irregular B-spline grid of 17*31=527 points is used for the background. The velocity of the salt body is constant (1 B-spline gridpoint). Total number of data points is thus 128+21412=21540 points, total number of B-spline points is 56+527+1=584 points. Program 'gmdmgmt.for' then requires 584*21540+21540+584*584=12941956 storage locations in array RAM and we may perform the calculations with MRAM=14000000, which is a value suitable for a computer with 64 Mbytes of memory.
To estimate the maximum sum of travel times from the source
and receiver,
we compute average slowness 0.120611 E-3 s/feet.
If we consider
maximum sum of ray paths about 50000 feet, we arrive at the estimation,
that the maximum sum of travel times from the source and receiver is
limited by Tmax = 6 s.
According to the
equation (2) of the paper
'mar-inv.htm',
we choose the maximum
Lyapunov exponent lambda of the new model as
lambda .LE. 2/Tmax = 0.333 /s . | (1) |
According to the
equation (7) of the paper
'mar-inv.htm',
making the same approximation as in equation (10) of the paper,
we arrive at
||v|| .LE. sqrt(8/3) * ( 4 * (1+ln(2)) / Tmax )2 / average(v) , | (2) |
||v||max = 2.6 E-4 (s*feet)-1 . | (3) |
In this example, we try to apply the minimum reasonable smoothing of the interfaces. We thus do not prescribe maximum Sobolev norm of the interfaces.
Objective function is taken in the form of
y = (|x-x0|/(1 foot))2 + (SOBMULx*||x||)2 + (|v-v0|/(1 foot/s))2 + (SOBMULv*||v||)2 , | (4) |
We use unit given velocity deviation and unit given interface deviation in the selected objective function. The unit deviations are realized by proper selections of parameter ERRMUL, which should be chosen as a square root of the ratio of the whole model volume to the volume corresponding to a single data point of the considered data set.
The inversion of the interfaces and the inversion of the velocity field are mutually independent in this example. We can thus perform the inversion of the interfaces in the first step, and the inversion of the velocity in the second step.
(1) Inversion of the interfaces:
We choose SOBMULv=0 s*feet, and perform the inversion without smoothing of the velocity, which is not of our interest at the moment.
In this example, we try to apply the minimum reasonable smoothing of the interfaces. We thus start the inversion without smoothing of the interfaces, with SOBMULx=0 feet (the interfaces are smoothed just by the projection onto the B-splines). We obtain oscillations of both surfaces, which produces other salt body than that one given by data. In the next steps of inversion we thus include minimization of the curvature of the interfaces by increasing the value of SOBMULx. When we reach approximately the value of SOBMULx=2500 feet, the spurious salt body moves outside the model volume. We thus consider the value of SOBMULx=2500 feet as the value corresponding to the minimum reasonable smoothing of the interfaces. Whether the interfaces are sufficiently smooth will be tested by ray tracing.
We compute standard interfaces deviation
|x-x0|=5.6 feet,
stored in the file 'hes-sd1.out'. This deviation causes
a travel time error dT approximately
dT = |x-x0| * |1/average(v) - 1/salt(v)| . | (5) |
(2) Inversion of the velocity field:
We have selected the value of SOBMULx=2500 feet, which ensures proper position of the interfaces. We will keep the value during the velocity inversion.
(a) First iteration:
We first perform inversion with SOBMULv=0 s*feet
(i.e., without smoothing),
and obtain the standard velocity deviation of the model
|v-v0| = 13.94 feet/s | (6) |
||v|| = 37 E-4 (s*feet)-1 , | (7) |
(b) Second iteration:
We choose
SOBMULv = |v-v0| / ||v||max , | (8) |
SOBMULv = 13.94 / 2.6 E-4 s*feet = 53615 s*feet | (9) |
SOBMULv = 50000 s*feet | (10) |
||v|| = 2.9 E-4 (s*feet)-1 , | (11) |
(c) Subsequent iterations:
We choose several values of SOBMULv for subsequent
iterations in order
to see the behavior of the model in dependance on SOBMULv.
After each inversion we perform several tests, which are described
below in the chapter
"Notes on file 'hes-test.h'",
and compute for each model
the value of average Lyapunov exponent lambdaav and
the maximum number of arrivals NUMmax from the source
located in the left top corner of the model.
We performed the following tests of SOBMULv:
model | |v-v0| / (feet/s) | ||v|| / ((s*feet)-1) | lambdaav / (s-1) | NUMmax |
desired values | 2.6 E-4 | 0.333 | 7 | |
model 'hess.dat' | 0.423 | 11 | ||
SOBMULv=0 | 13.94 | 37 E-4 | 0.380 | 7 |
SOBMULv=17500 | 15.19 | 4.9 E-4 | 0.345 | 6 |
SOBMULv=35000 | 17.41 | 3.5 E-4 | 0.316 | 5 |
SOBMULv=50000 | 19.24 | 2.9 E-4 | 0.300 | 5 |
SOBMULv=70000 | 21.42 | 2.4 E-4 | 0.287 | 5 |
We can see, that for SOBMULv=0 s*feet we obtain the standard velocity deviation about 14 feet/s. Desired value of lambdaav=0.333 s-1 corresponds to the velocity deviation about 16 feet/s. Most of the velocity deviation is thus caused by the projection onto the B-splines, and we should try the inversion with finer parametrization of the model.
A figure of Fourier spectrum of velocity deviation should help us
to decide, in which direction and how many times we should
densify the parametrization of the model. We must keep in mind,
that the parametrization must not be denser than the data.
Fig. 2: Fourier spectrum of velocity deviation in the model obtained with SOBMULv=50000 s*feet. |
We choose new parametrization of the model as follows: for the interfaces we keep the irregular B-spline grid of 28 points; for the background velocity we take horizontally regular grid of 26 points with spacing 1600 feet and vertically irregular grid of 58 points with spacing from 200 to 400 feet. The parametrization may be found in file 'hes2-mod.dat', which is used as a new starting model for the inversion.
We have now 56+1+26*58=1565 B-spline points, and the memory requirements are 1565*21540+21540+1565*1565=36180865 storage locations. This is too much and we thus cannot use all the data. We resample the data in horizontal direction and use the data grid of 67 points with spacing 600 feet horizontally and original 106 points with spacing 200 feet vertically. We have now 67*106+2*64=7230 data points and the memory requirements are 1565*7230+7230+1565*1565=13771405 storage locations.
We performed the following tests of SOBMULv with new starting
model 'hes2-mod.dat':
model | |v-v0| / (feet/s) | ||v|| / ((s*feet)-1) | lambdaav / (s-1) | NUMmax |
desired values | 2.6 E-4 | 0.333 | 7 | |
model 'hess.dat' | 0.423 | 11 | ||
SOBMULv=1000 | 5.52 | 12.3 E-4 | 0.398 | 6 |
SOBMULv=10000 | 7.28 | 6.5 E-4 | 0.365 | 5 |
SOBMULv=17500 | 9.03 | 5.2 E-4 | 0.343 | 5 |
SOBMULv=25000 | 10.68 | 4.4 E-4 | 0.327 | 5 |
SOBMULv=35000 | 12.75 | 3.7 E-4 | 0.313 | 5 |
Note that with the new parametrization we are not able to perform the inversion with SOBMULv=0 s*feet, because we have now more background velocity B-spline gridpoints in the salt body, and the inversion becomes ill-conditioned.
History file 'hes-test.h' serves for numerical tests of model 'hes-mod.out' obtained after the inversion realized by history file 'hes-inv.h'. The figures displayed below are obtained by running 'hes-inv.h' with starting model 'hes2-mod.dat' and the values of SOBMULx=2500 feet and SOBMULv=25000 s*feet.
We discretize the velocity
in the smoothed model 'hes-mod.out'. We calculate
the gridded velocity deviation from the velocity in the model
'hess.dat', the gridded relative
velocity deviation, the
standard velocity deviation and the standard relative
velocity deviation, the Fourier spectrum
of the velocity deviation, the deviations of the interfaces,
the relative deviations of the interfaces, and the standard
deviation of the
interfaces
and standard
relative deviation of the interfaces.
Fig. 3: Velocity in the smoothed model 'hes-mod.out'. |
Fig. 4: Velocity deviation. |
Fig. 5: Relative velocity deviation. |
Fig. 6: Fourier spectrum of velocity deviation. |
Directional Lyapunov exponents and the
average Lyapunov exponent for the 2-D model are calculated by program
'modle2d.for'.
Fig. 7: Angular dependence of directional Lyapunov exponents and the average Lyapunov exponent. |
Three tests of ray tracing of a refracted wave
followed by interpolation between rays are performed using programs
'crt.for' and
'mtt.for'. The three sources are located
at the top of the model in the corners and in the middle of the model.
Fig. 8: Rays and numbers of arrivals for the source in the left top corner of the model. |
Fig. 9: Rays and numbers of arrivals for the source in the middle of the top of the model. |
Fig. 10: Rays and numbers of arrivals for the source in the right top corner of the model. |