OCTAVE CAINV3

Down-hole SH-wave data are inverted for stiffness and damping with this Octave program. Required are the *.his file results from programs BVAS (section 8.2.5) and BAMP (section 8.2.6). Files bvas.his and bamp.his are required for each depth interval of interest. In order to compute uncertainty error bars, the depth interval should include as many subsurface stations as possible. Since cainv3.m is a joint inversion of body wave velocity dispersion and amplitude decay (corrected for beam divergence), the *.his files do not need to include exactly the same subsurface stations, as when there is a need to remove poor data from one or both *.his files. A companion program that computes the forward problem is cafwd3.m see section 9.1.1 .

The governing differential equation for this problem is a 3rd order PDE that is formulated as a 1-D plane wave problem (hence the need for the BAMP program to correct for beam divergence.

$$\frac{\partial^{2}u}{\partial t^{2}}=C_{1}\frac{\partial^{2}u}{\partial x^{2}}+C_{2} \frac{\partial^{3}u}{\partial t\partial x^{2}}$$ (4)

where “u” is particle displacement, “t” is time, “x” is the coordinate in the direction of wave propagation, $C_{1}$is the stiffness $\left(\frac{m^{2}}{s^{2}}\right)$, and $C_{2}$ is the damping $\left(\frac{m^{2}}{s}\right)$. Equation (4) reduces to the elastic wave equation when the damping value, $C_{2}=0$. In that case, the phase velocity is constant for all frequencies, and the wave does not experience any decay (for a 1-D plane wave). In the elastic case, the phase velocity will be $\sqrt{C_{1}}$ .

In the more general case, $C_{2}\neq0$, and there will be both velocity dispersion and exponential, inelastic amplitude decay. A solution of equation (4) is

$$u(x,t)=\exp\left(-\alpha x\right)\cdot\cos\left(\beta x-\omega t\right),$$

where the wavenumber is complex and given by $\beta+i\alpha$ . Michaels Michaels (1998) shows that the inelastic decay of a plane wave will be given by

$$\alpha=\frac{4\sqrt{D}\omega^{2}C_{2}}{\left(2\omega C_{2}\right)^{2}+D^{2}}$$

where $\omega$ is angular frequency (rad/s) and the quantity, D, is given by

$$D=2\left(C_{1}+\sqrt{C_{1}^{2}+\omega^{2}C_{2}^{2}}\right).$$ (5)

The phase velocity, c, varies with frequency according to the following relationship

$$c=\frac{2\omega^{2}C_{2}}{D\alpha}.$$ (6)

The values for $C_{1}$ and $C_{2}$ can be expressed in terms of the following :

$$C_{1}=\frac{\left(\beta^{2}-\alpha^{2}\right)\omega^{2}}{\left(\beta^{2}+\alpha^{2}\right)^{2}},$$ (7)

and

$$C_{2}=\frac{2\alpha\beta\omega}{\left(\beta^{2}+\alpha^{2}\right)^{2}}.$$ (8)

Determination of $C_{1}$ and $C_{2}$ is by nonlinear joint inversion of the phase velocity, c, and inelastic decay, $\alpha$ , over a range of frequencies. The inversion is currently performed in the Octave procedure, cainv3.m. Initial estimates of stiffness and damping are obtained at the frequency corresponding to the largest $\alpha$ measured by bamp. First, $C_{1}$ is found by evaluation of equation (7). In that computation, $\beta=\frac{\omega}{c}$ . Then, $C_{2}$ is estimated from equation (8).

RUNNING CAINV3:
Start an octave session, type cainv3
GUI, use mouse to pick min and max frequencies, click OK
and then use the mouse.  Horizontal position is all that
is read.  Focus one of the panels.
You can exclude some frequencies, and that will create
an fbx vector.  If you include all frequencies, you may get
an error statement (since it can't write out something that
does not exist).  Typically not a problem when you run
caplot3.m later.  Don't worry about it.

GUI, C1=stiffness, C2=damping initial estimate for the
3rd order wave equation.
GUI, Choose weighting
GUI, Choose balance between damping and velocity, .5 good idea
Plots, update as inversion progresses
GUI, continue LSQE plot
GUI, continue Chi squared plot
GUI, save results to disk, yes if you want to run caplot3.m

Figure 31: CAINV3: First display. Use mouse to pick frequency limits for analysis, low and then high.

After running cainv3, you may wish to make nice plots. For this, there is program caplot3 (8.2.8).