More Visualization with Mathematica

Introduction and Background

As we learned earlier, visualization is a key ingredient in the
design and analysis of engineering systems.

Throughout history, scientists and engineers have been using visualization
as a way to better understand their models and analyze their solutions.  In
the beginning, and for some time, engineers used hand drawn plots to
graphically illustrate their problems.  While this worked fairly well for
static, two-dimensional phenomena, it was often very difficult---if not
impossible---to graphically illustrate dynamic and/or three-dimensional
phenomena.  Often,  scientists and engineers would have to resort to
actually building a small scale model of the engineering system before they
could grasp some of the inherent complexity and begin to understand the
system.

The next step to visualization came about when large mechanical plotters
were developed.  These plotters could graph complex two-dimensional data
and allowed scientists to ``slice'' through three-dimensional data.  While
engineers could now obtain some clues into understanding their complex
data, it wasn't until the advent of three-dimensional computer graphics
systems that scientists and engineers could really begin to fully analyze
complex three-dimensional and/or dynamic systems.  It may interest you to
know that some of the world's very first visualizations of scientific and
engineering data occurred here at the University of Utah.  For example, it
was at Utah that the visualization of complex molecules and results from
finite element analysis were first realized.

When one thinks about scientific visualization, one usually thinks about
ways to view the geometrical model and/or interpret and represent the
computer solution which approximates some real world experiment or design.
And while that is our primary goal in this lesson, we should note that
visualization can also be used in other aspects of engineering analysis,
including checking data entry and integrity.  In past two decades,
computing power has increased at an exponential rate such that we can now
simulate highly complex phenomena using millions of pieces of information.
Visualization can be used to check the integrity of our input (and output)
data.  Furthermore, graphical means can also be used to enter data as it is
not often feasible (or desirable) to enter large data sets by hand.
Visualization can also be used in the analysis of program performance and
aiding to debug Fortran and/or C code.  To analyze large, complex system it
often requires a large number of instructions in the form of Fortran/C
code.  Visualization systems have been (and are being) developed which
allow engineers to find places in their code that are inefficient and/or
are in error.

Earlier we investigated how to use some of Mathematica's two-dimensional
plotting feature.  In this lesson, we will experiment with more advanced
two-dimensional plotting as well as three-dimensional plotting of an
animation.

Population Growth: Graphing Tabulated Data

One subject that crosses many disciplines is that of modeling population
growth and or decay.  Biologists are interested in such models in terms of
modeling species interactions, economists are interested in such models for
use in economic forecasting (as we saw in a previous lesson on compound
interest), bio- and chemical-engineers are interested in such models for
predicting the rates of certain reactions.  Furthermore, the same
mathematical model can be used to describe radioactive decay, diffusion of
a chemical across a thin membrane, the rate of cooling of objects, and even
the spread of disease.

Let's look at a simple population growth model.  Let P(t) be the size of
the population at a given time t.  Let b and d be the birth and death
rates per unit time, respectively.  It seems reasonable that the change in
the population over time would be proportional to the number of births
minus the number of deaths over time.  Or

DP/Dt is proportional to  (b-d) Dt

where

D

signifies a discrete change.  If one sets these two
expressions equal and takes the limit as
Dt goes to 0 (remember
your basic calculus), we can write this relationship mathematically as,

dP(t)/dt = (b-d) t

If one integrates this equation, one can find a solution for the population
at a given time as a function of the number of births and deaths,

P(t) = P(0) e^((b-d)t)

where P(0) denotes the population at time equal to zero (i.e.,  whenever
you start measuring it).

As we can easily see, the population will grow at an exponential rate if
b-d > 0, decay at an exponential rate if b-d < 0, or stagnate if
b-d = 0.

Let's assume we have a small town which, because of its scenic beauty and
appeal to rugged outdoors types, has been experiencing steady growth over the
last several years.  We can plot the town's growth using Mathematica's ability to
graph sets of tabulated data.

Using emacs, create a file called data.m containing the following data
points recorded from census data over the last ten years.  The numbers in the
left column signify the year and the number in the right column signify the
population at a given year.

The data was created in Mathematica by evaluating equation (7.3) with an
initial population of P(0) = 1000 and a birth minus death constant of
.2 .

Now you can read this information into Mathematica with

<<data.m

The ListPlot function will plot from tabulated data in this
form.  Try

ListPlot[data];

You will notice that Mathematica does not  ``connect the dots,''
for you.   Whenever you plot a set of tabulated data,
Mathematica automatically defaults to this value. Instead we could
have Mathematica find the ``best possible'' smooth curve through the data 
points and plot that.  How one finds this ``best
possible'' curve through the data points is an important topic which
we will revisit in another lesson.  For now we will just connect the
data points using straight lines between each point.  Try the following:

ListPlot[data, PlotJoined-> True];

As you might guess, he PlotJoined option tells the Mathematica 
ListPlot function whether or not the data points should be joined with a 
straight line.

Animation and 3D Plots

Beyond producing static two-dimensional plots, Mathematica can produce three
dimensional plots and can animate series of two and
three-dimensional graphs.

Briefly, most of the animation features can be utilized directly in the Mathematica
``front end''.  Additional features can be found in a set of  functions
contained within the Graphics/Animation package library.  All of these features
can be used in conjunction with the various plotting functions

Mathematica comes with a very  large set of additional functions which 
are included in the form of ``packages'' that you load as you need them.  
The set of packages includes those dedidcated to 
Linear Algebra, various aspects of graphics, Statistics,
and several more.  By invoking the
command,

Needs["Graphics`Animation`"]

we can utilize the animation (and other)
functions contained in the Graphics/Animation library.
Note that once we have invoked this command
it is valid for the rest of our Mathematica session.  To find out more about
the relevant commands, use either the Mathematica book or 
the on-line help facilities of Mathematica to 
learn about the Needs function (For example, but inputting

?Needs

in a Mathematica session.

Also, if you are using a notebook interface to 
Mathematica you can see what functions are defined

To try out the animation and 3D graphics features, let's consider another
physical phenomena which crosses many disciplines, that of wave propagation.
In general there are two types of waves,  transverse and 
compressional.  A  transverse wave occurs when the impulse starting the
wave propagation is perpendicular to the direction of the resulting traveling
wave.  Examples of transverse waves include, waves on a string (like plucking a
guitar), water waves, and light.  A  compressional wave occurs when the
impulse starting the wave is parallel to the direction the wave travels.
Examples of compressional waves include sound and earthquakes.

Let's consider the motion of transverse waves.  For those of you that are
interested, the equation for transverse waves is given as,

where u is the position of the wave as a function of space and time and
v is the velocity of the wave.  While the equation might look complicated
at first, it actually has a fairly simple solution which can be expressed
in terms of a cosine.  For a one-dimensional wave we have,

u(x,t) = A cos((2 p/l) (x - vt))

where A is the amplitude of the wave (i.e., the ``height'') and
l is the wavelength (the length of a wave is measured from crest to crest---or
from any two like points).  Depending on our choice of the amplitude and
wavelength, the solution should look very familiar to us---like a cosine wave
that moves in time.  We can utilize Mathematica's animation features to see this.
Let's set A=1.0,l=5.0, and v=2.0.  Then the wave is a function of two variables, 
x and t, given by

u[x_,t_]:=Cos[(2 Pi/5.0)*(x-2.0 t)];

Now we make a ``table'' of plots of this function using the Animate function 
which we loaded via the Graphics/Animation package

Animate[
		Plot[u[x,t], {x, 0, 10}], 
		{t, 0, 2Pi, Pi/8}
	];

This is the solution of a wave traveling along a string (the type of
wave you would create by taking a rope and moving it up and down or by
dropping a stone into a calm pond).  If you are using the Mathematica front 
end for this lesson you should be able to select the collection of plots and
make them ``scroll'' on the screen like a cartoon. If you run this in a constant
loop, it may remind you of the output of an oscilloscope (or a scene
from an old science fiction movie) - depending on our fast your
computer is processing the information.

Now experiment with different values for A, v, and
l to see how they influence the shape of the wave. It will be more 
convenient if you make a new function, wave, that explicitly makes all
of these parameters variables.  For example,

wave[x_,t_,a_,v_,lam_]:= 
a Cos[(2Pi/lam)*(x-v*t)];

While the solution for the traveling wave along a string isn't too complex,
things change when the waves are moving in two and three dimensions,
often times with different amplitudes and wave speeds along each direction.
Let's look at a couple of waves at different instances in time for a
two-dimensional wave. Let's start off with the simplest case of a wave in the
x direction and one in the y direction which have exactly the same parameters:

Plot3D[
	wave[x,1,1,2,5]+wave[y ,1,1,2,5], 
	{x, 0, 2Pi}, {y, 0, 2Pi}
];

This is the type of wave motion you would obtain by moving a piece of
cloth---such as a sheet---up and down).  For the above parameter values,
the wave should appear symmetric. Notice how difficult it would be to 
reconstruct such a surface just from a set of numbers or two 
dimensional plots.  

Usually it is important that the values for v and l in the two functions above match.  
The velocity and wavelength are physical parameters related to the 
material in which the wave is propagating.   If these parameters are not the same 
in each direction then usually that means you have made a mistake, since otherwise 
what we are plotting would not make sense physically. However, there are some 
substances, called anisotropic materials, which do have different wave velocities
and wavelengths in different directions.  We will not deal with these here, so for us
v and l should always agree. 

Using this same physical problem, we can
illustrate another type of visualization. Mathematica can also make contour 
plots of two dimensional data.

ContourPlot[
	wave[x,1,1,2,5]+wave[y ,1,1,2,5], 
	{x, 0, 2Pi}, {y, 0, 2Pi}
	];

(notice that as our Mathematica commands get longer and more complicated
we use tricks like indenting to make what we are doing a little clearer)

What you should see is something that reminds you of a topographical
map.  The contours represent lines of equal values of the function, i.e., places of
the same amplitude.  This is a direct analogy to a topographical map
which represents lines of equal elevation.  Similar graphs can be used
to show places of equal temperature, or equal voltage, etc.
Iso-contour maps are used in almost every field of science and
engineering.

Let's look at another time instance; notice how the pattern has shifted.

ContourPlot[
	wave[x,2,1,2,5]+wave[y,2,1,2,5], 
	{x, 0, 2Pi}, {y, 0, 2Pi}
	];

Comparing this to our previous contour plot an interesting pattern should
emerge.  Just as strings have different frequencies which give a
so-called  standing wave at which they vibrate (these are called
 normal modes).  Vibrating surfaces have many more different
combinations of frequencies at which they can vibrate and produce what
are called nodal lines.  Along these lines there is no motion of
the surface.  An 18th century scientist, Ernst Chladni, figured this
out by placing sand on the surface of a flat plate and then vibrating
the plate with a violin bow.  This caused the sand to distribute along
the nodal lines and form interesting patterns.  What you are seeing is
the two dimensional projection of these patterns.   Note: If you
want to try this at home, just get a flat metal plate and some sand
and vibrate the edges with a bow or tightly held cord.

By using three-dimensional graphics, we can look at the nodal lines as they
appear throughout the surface and also, we can animate the wave propagation in
three-dimensions.  Let's try a  three-dimensional animation!  We substitute
the variable t into the above function, and again make a ``table'' of plots for 
different values of t.

Animate[
	Plot3D[
		wave[x,t,1,2,5]+wave[y,t,1,2,5], 
		{x, 0, 2Pi}, {y, 0, 2Pi}
		],
	{t, 0, 10, 2}
	];

Where we have made t a variable which we allow
to vary from 1 to 10.  Now experiment with the values
for l in the original 3D plot and see how the waves (and nodal
patterns) change.

Continue to experiment and change the amplitudes and velocities of the wave
and notice their effects.  You are now practicing a form of computational
science---you are using scientific visualization to investigate complex
wave phenomena that would be otherwise impossible with pencil and paper or
sets of data points.  Hopefully this will give you a feeling for the power
of scientific visualization and some of its many uses.

Homework

Problem 1: In the population growth problem, if 
b-d > 0, the population will double its present size at time T,

2P(0) = P(0)e^((b-d)T)
 
 or
 
 e^((b-d)T)=2

or

T = ln(2)/(b-d)

where T is called the doubling period.  If the world currently has a
population of 4.5 billion people and its doubling period is 30 years
(i.e., T=30), graph the population of the world after 30, 50, 100
and 200 years from now on the same graph.  Hint: Utilize equation
(7.3).

Just for fun (not for credit), the surface area of the earth is about
1,860,000 billion square feet, how much space would each person have after
1000 years?

Problem 2: In the examples of wave motion that we studied in the
lesson, the waves traveled at a constant velocity.  Such waves are
called nondispersive due to the linear relationship between the
speed of the wave and its length.  This means that they don't ever
slow down, but always maintain the same speed.  In the real world, the
velocity of the waves changes as a function of the properties of the
media that it travels through, primarily due to friction.  Such waves
are called---you guessed it--- dispersive waves and the form of
the wave changes over time.

(a) Given the original description, the velocity of the wave, v, was set as
a constant value.  Now we will modify this to allow the velocity to be a
function of the wavelength, where

v(2p/l) = a + b(2p/l)

Let  a = 1.7 and b = (0.2).  Redefine the equation for the
one-dimensional wave u such that the new value of the velocity has
been incorporated.  Hint: You can use the definition that we already
gave for the one-dimensional wave and modify it. You can use the same
values for A and l as in the example.

(b) Now construct a graph which has the original expression for the
one-dimensional wave that was illustrated in the lesson  and also
includes the modified wave expression from part (a) on the same plot,
for a value of t = 1.0  Are there any differences?  What about at
t=3.0?

(c) In order to investigate any possible difference further, 
animate the plot you made in (b) of both waves as a function of time
on the same graph.  Use t = 1 up to 20 or 30 and approximately 50
frames.   Describe what happens and how.

(d) Now repeat steps (a) and (b) for the two-dimensional wave and
report what you observe.  Hint: Again, you can build directly from the
examples given in the lesson.

(Just for fun) You can compute the phase velocity of the wave by
graphically comparing the distance the wave travels from one ``crest'' to
another.  You can also measure the ``dispersiveness'' by comparing the same
instance of time between the dispersive wave and nondispersive wave on the
same graph.

Prepare a report giving the answers to these questions.  
You don't have to hand in any plots, but be sure to tell us
how the plots and animations that you did helped you to answer the questions.
Submit your report for grading.