Geol 455/655 - Lab 2, Courant Condition
For the second exercise, the objective is to begin exploring
the parameters of the numerical simulation.
The first one is the Courant Condition, which affects the
Aliasing and the Courant Condition
The blue curve above is a sine wave with wavelength lambda.
On the left, the green lines represent a discrete
approximation to the continuous curve, with two sample
points (white circles) per wavelength. The green lines
are clearly an adequate approximation of the sine wave,
with insignificant differences between the two.
On the right the red lines represent a sampling of the
sine wave at only 1.5 points per wavelength. This is
clearly not enough samples per wavelength, as the
red lines do not approximate the blue curve well at all,
and there are large errors between them.
In E3D's computations, the time-step dt must be short
enough to properly represent the waves being modeled.
The equation v = f (lambda) above shows that the
velocity at any point in the model is related to the
The Nyquist Criterion states that
a discrete representation of a wave in time must have
a sampling frequency fs
that provides two points per period, or
fs = 2
fmax , where
fmax is the highest
frequency found in the waves,
fs = 1/dt ,
and dt is the time sampling interval.
The Courant Condition states that the time sampling
dt must be small enough that the longest wavelengths,
propagating at the highest velocity
Vmax , do not outrun
the spatial grid sampling dh. There is an additional
adjustment to the Courant Condition that we will not
worry about here, since it comes from the nature of
the finite-difference numerical approximation used by E3D.
Combining the Nyquist Criterion with v = f (lambda)
will show you how small you have to make dt to keep
two samples per spatial wavelength.
- Re-run the computation from Lab 1, to make sure you
still have all the input and output files in your
folder. As you watch it run, note what E3D says about
the Courant Condition, ct/dt, and the MFLOPS achieved.
- Re-name the sac files from the Lab 1 computation,
to indicate the value of the dt parameter used on the
``time'' line of the ``test.in'' and ``test1.in'' files:
cp LSM.000.x LSM.12dt.x
cp UNLV.001.x UNLV.12dt.x
- Edit the ``test1.in'' file to change the ``dt'' and ``t''
parameters, making dt half of what it was, and doubling
t so the same amount of model time is run. Don't bother with
the ``test.in'' file, because for this lab we will not
need to run ModelAssembler. The time line will look like:
# time dt=0.12 t=1000
time dt=0.06 t=2000
So the previous values are just commented out. This should
provide finer time sampling, and lead to more accuracy in
- Run the model with the command:
/quake/rg/e3d/bin/e3d test1.in &
We don't need to run "csh run.sh" because we already
have all the medium volume files, and don't need to run
- Observe the wave propagation in the pop-up image window. What is different about this run compared to the Lab 1
run? The Courant Condition? How much longer (wall-clock time)
does it take? Why doesn't the MFLOPS change?
- Re-name the sac files from this new computation,
to indicate the new value of the dt parameter:
cp LSM.000.x LSM.06dt.x
cp UNLV.001.x UNLV.06dt.x
- Use ``sac'' to plot and compare the four seismograms:
the two LSM and UNLV x-components from the Lab 1 run
with dt=0.12 s; and the two with dt=0.06 s. In sac you
can get them with the command:
read LSM.12dt.x LSM.06dt.x UNLV.12dt.x UNLV.06dt.x
Note that the seismogram files made with dt=0.06 s take up
twice as much space as those made with dt=0.12 s, and
they took longer to compute. Since they all cover the same
120-sec time period, sac plots all four seismograms on the
same horizontal axis.
Are the more expensive seismograms using dt=0.06 s any
different or better?
- Now, just as in preparing for the dt=0.06 s run as
above, edit test1.in and make another E3D run using
dt=0.24 and t=500. What happens? What is the Courant
Condition for this run?
What to turn in:
- Answer: Given constant frequency, how do wavelengths
change as velocity changes from place to place in the
grid? Using this concept, explain how wavelengths on the
surface of the computation change across the visual display
while E3D is running.
- Answer: Given the usual relation between P and S velocities,
will P waves or S waves have longer wavelengths?
- Answer: Solve v = f (lambda) for
fmax given the
maximum velocity reported in the grids by E3D, and
using dh for lambda. According to the Nyquist
criterion, what dt is required? Is this similar to E3D's
reported dt as required by the Courant Condition?
- Answer the questions in italics above.
- Turn in the plot of the 4 seismograms from the
2 stations and the 2 different dt values, properly and
completely labeled on all axes (with a pencil if necessary).
- Answer: Given velocity grids and a fixed dh such as we
are using here, could we save any computational effort by
removing structural features such as the basins? Can you
propose a feature of the models that we could remove to
enable us to make dt larger, and thus make a synthetic
seismogram covering the same time with less computational
effort? (Are there any structural features in the grids
other than the basins?)