CBS 613: COMPUTER LAB I
In future labs we will be using different software to compute numerical approximations of the solutions to ordinary and partial differential equation models of certain biological systems. This lab is meant to provide an overview of some of the software available and what kinds of things they can do.
The math department uses a UNIX operating system for its computers. UNIX is very powerful but hard to use at times. For a (very) brief introduction read thru the following link:
Changing you password
Think of a new password for your account. Make it 8 characters. It must include a non-alpha character or else the system will reject it. Then in the terminal window type:
passwd
You will be asked to type in you old password and your new password twice. Then next time you log in,.....remember the new password!
Some Software:
We will use (at most) three software programs called "matlab", "xppaut" and "xtc". All can (numerically) solve Ordinary Differential Equations (ODE) and Partial Differential Equations (PDE). To run both "xppaut" and "xtc", a text file must first be created which defines the problem, parameter values and the like. For the most part, "matlab" does not require such a file but it can read special files whose suffix is "*.m" (as in file.m, bob.m, problem.m). This becomes convenient when your problems get more complicated. Both "xppaut" and "xtc" have their own suffixes.
| Program Name | Suffix | UNIX command to execute |
| matlab | *.m | matlab -nodesktop |
| xppaut | *.ode | xppaut file.ode |
| xtc | *.xtc | xtc file.xtc |
The mathematics Department has different versions of matlab: Version 5.3 and Version 6.1. Although they use the same commands, Version 6.1 has a fancy window interface. I find this version buggy and somewhat confusing to use. So, for the purpose of this lab you will introduced to the older version 5.3. To access this older version, in your UNIX terminal window, type (case sensitive):
alias Matlab /opt/matlab5.3/bin/matlab
Matlab
and wait for the program to load. You should have a new ">>" prompt. Any command entered after this prompt must now be a matlab command.
Type in a few simple commands:
>> x=0:10
>> y=(x^2)/10
>> plot(x,y)
Yah,....so what? My calculator can do that! So what makes matlab so useful is that it has a huge number of commands for doing virtually every type of mathematics and statistics your average scientist will ever need. To get an idea of how big matlab is, let's spawn the helpwindow....type:
>> helpwin
You should have windows that look like:
As you scroll down this list, each item has sublists (double-click on them) of items related to calculations, mathematical algorithms, graphing, programming computer interfaces.... From this help window you can explore and learn. Just to get a taste, lets try some demos:
I've picked out a few in the following table. Some of the demos are informative. Others are "cute" and not so informative. The thing to remember is that everything you see is created by a sequence of commands you can type in the ">>" prompt or enter into a *.m file. Do them in sequence, read them and try to observe the code that results in what you see. (I will interject to make sure you understand)
| Demo Name | Matlab Command to Execute | Description |
| intro | intro | basic matrix and vector calculations |
| fitdemo | fitdemo | nonlinear curev fitting example |
| census | census | curve fitting demo for US census data |
| funfuns | funfuns | minima, roots, integration |
| lotkademo | lotkademo | Demo which solves Lotka-Volterra ODE |
| graf2d2 | graf2d2 | Demonstration of some 3-D graphics |
| xfourier | xfourier | Demonstration of Fourier Series Approxmations |
| odedemo | odedemo | Multiple ODE models and methods demo |
Solving planar ODE's with "pplane" "pplane" is another matlab command which is useful for solving systems of ODE's with two dependent variables (planar systems). We will demonstrate its use on a model of a set of chemical reactions.
Consider four chemical species X,Y,A and B with the associated chemical kinetics:
Since X is involved in its own production, the reactions are autocatalytic. Using the Law of Mass Action , differential equations for the (nondimensionalized) concentration x of X, y of Y, a of A and b of B can be derived:
Here t is time and the "concentrations" are scaled with respect to combinations of the rate coefficients ki, i=1,2,3. For the purpose of the following simulations we will regard the initial concentrations of A and B being fixed. This would be a reasonable approximation if there were many more A and B molecules than X and Y, for instance.
In the matlab window, type:
>> pplane
Alter (by clicking/backspace) the differential equations and parameter values so that you get the following:
When you are done, click "proceed". You should get:
The arrows indicate the vector field for the differential equations. To solve the differential equations you need to stipulate initial conditions , i.e., the initial concentrations x(0) and y(0). To do this, point the mouse arrow at a point and then click. Do this a bunch of times. You should get something like
For these a and b, the concentrations x(t) and y(t) approach constant values as t increases (approaches an equilibria).
Go bakc to the other window and change a and b to:
a=0.1
b=0.5
then click "proceed" again. Try a few new initial conditions and you should get oscillations.
If you do not want to view the concentrations as a phase portrait but would rather see x and y versus time,
1) Select "Both" under the "Graph" menu
2) Point the cross hairs at one of your orbits and click on it
You should get a new window that looks something like :
Go back to the phase portrait window, then:
3) select "Erase all solutions" under the "Edit" menu
4) select "Show nullcines" under the "solutions" menu
5) select "find an equilibria" under the "solutions" menu
6) point the cross hairs at the intersection of these nullclines and then click.
You should get a new window:
The window above prints the results of a (numerically computed) linear analysis of the equilibria. The result states the equilibria is an spiral source - hence unstable.
pplane and matlab does lots, but we'll stop for now. To exit matlab, type:
>> quit
in the original matlab window.
Unlike matlab, "xppaut" is specialized to solve and analyze the dynamics of ODE's and Maps. The program was written by Prof. B Ermentrout (U. Pitt, Math Dept) and, unlike matlab, it is free. It comes with a full web based tutorial and PDF documentation. As mentioned earlier, the user must first create a file with a *.ode suffix which defines the problem. Click on the link below. The code defines a set of differential equations for an model of electrical activity (transmembrane potential) in cells whose channel dynamics are slowly modulated by calcium. I will briefly describe the code before you save it:
Burst.ode Bursting Oscialltor
Before, let's take a quick look at the web based tutorial:
We won't examine this program in detail now. For today, we'll just do a couple of simple things.
In the UNIX terminal window, enter
xppaut Burst.ode
Tons of windows come up. Close all but the black window. Then
A) select "Go" from the "Initial Conditions" menu.
You should get something like
This is a numerical solution the model equations (voltage u(t)) for the initial conditions, numerical parameters and model parameters defined in the Burst.ode file. These you can change interactively using the menu items. Try this:
B) click on "Parameters" in the black window.
C) at the top of the window type "eta" then hit the return key
D) Set eta=1.5, then "return"
E) Click "Erase"
F) Select "Go" from the "Initial Conditions" menu
Note in this case how the qualitative nature of the solution changes. To exit:
Select "Quit" under the "File" menu.
We'll do this at a later date