Computer Science Theory tells us:
Repetitive operations "R" computers. E.g.:
Counted loops
(
primitive recursion in recursive function theory):
perform loop some given number of times
(e.g., "brush hair 100 strokes before bed.").
Conditional loops
(
general recursion):
repeat until done
(e.g., "whip until fluffy")
In Matlab, since matrices are the primary datatype, explicit looping is often not needed. Vectorized commands operate on entire matrices. With Matlab we take this for granted, but it's both a luxury and a curse: a powerful capability that can be misused.
Attaway often refers to vectorized versions of commands as
"efficient". They are terse to us humans, but in fact
vectorized commands are implemented under the hood as loops,
so they can't be inherently more efficient in terms of execution time.
For example,
avec = sqrt(bvec*5 + 4);
is_greater_vec = some_vec > 5.0;
% produces logical vector
Good Rule: create a matrix or vector if you really need all its values (say to make a plot or as input for further processing). Do not construct a vector to hold input that can be generated and used one piece at a time.
For instance, to add all the numbers from one to a billion, we don't need a billion-long vector[1,2,3...]. We need one variable to hold the sum and we need to count: to generate the numbers to be summed one at a time. We need a for loop.
The vectorized
sum(1: 1000000000)
(sum one to a billion)
gives an out-of-space error, but
total = 0;
for value = 1:1000000000
total = total + value;
end
after running for several seconds, actually produces an (approximately) correct result.
Sometimes called a "do loop" (for Fortran), a
for loopvar = range_expr
statements;
end
Executes the statements (down to the range-expr
.
loopvar
.
range-expr
must be a range: one of the
simplest looks like 1:N, and the most trivial looks like
3. The trivial case only "does the loop once", so we don't
need
the loop!
for loopvar = x*2 + 10
statement;
statement;
end
executes the statements once with the value of the expression.
The range expression can be pretty general: Possibilities include operator expressions and vectors.
for loopvar = 1:N
for loopvar = 10:2:20
for loopvar = x: 2*y: z+150
for loopvar = [2 5 -98 13 pi]
% sets loopvar to each element
for loopvar = some_vec
% sets loopvar to each element
The index variable is so-called because it is so often used as an array index, in essence to find the address of the value we want.
i
and j
are generically used for index variables(i
and j
are
pre-defined constants that represent sqrt(-1). If you
assign to these variables, you will lose Matlab's intended semantics for them.
Enthusiastic hello:
for repetitions = 1:3 % i is set from 1 to 3
% but not used in loop!
disp('Howdy!');
end
Accumulator problem: accumulate partial answers into a single answer variable (very common).
For example, compute sum of n^{2}, as n goes from 1 to 10.
sumsq = 0; % a box for sum, initialized to 0
for n = 1:10
sumsq = sumsq + n^2; % update and save sum
end
What happens if we try the following?
for n = 1:10
sumsq = 0; % a box for sum, initialized to 0
sumsq = sumsq + n^2; % update and save sum
end
Above runs fine, just probably not what we want. So that's a problem... no error message.
On the other hand if we make up our own syntax and come out with gibberish like:
y = for x < 5
sum +x
end
We see an unhelpful
Illegal use of reserved keyword "for"
.
Which is technically right but
pretty vague. And since there are other problems in our little loop
above, either fixing one problem (or making it worse!)
is only going to lead to another useless error message. You're dying
here.
Moral: Trying to rewrite code without understanding what has gone wrong (aka Programming at Random) usually doesn't work.
More horrible examples of how NOT to Compute sum of n^{2}, as n goes from 1 to 10. Try to figure out what each of them will do...
sumsq = 0; % a box for the sum, initially 0
for n = 1:10
% n^2; % no
% sumsq = n^2; % no
% sumsq(n^2); % no! Stop Guessing!!!
end
Assume the following script is stored in for_script.m.
vecsum = 0; % a sum accumulator
vecprod = 1; % a product accum. (why not 0?)
N = length(a_vector);
for indexvar = 1:N % indexvar goes from 1 to N
vecsum = vecsum + a_vector(indexvar);
% add into accum.
vecprod = vecprod * a_vector(indexvar);
% mpy into accum.
end
vecsum
vecprod
Then
>> a_vector = 1:2:9; % = [ 1 3 5 7 9]
>> for_script
vecsum = 25
vecprod = 945
Several ways NOT to get the sum...
vecsum = 0;
N = length(a_vector);
for indexvar = 1:N % indexvar goes from 1 to N
% vecsum = sum; % syntax error
% vecsum = sum(1:N); % computes wrong sum N times
% indexvar; % what???!!!
% vecsum + indexvar; % computed and discarded.
% vecsum = vecsum + indexvar; % wrong sum
end
True, Matlab has built-ins for vector sum and product, but we're learning for-loops not memorizing a million Matlab commands.
Study the following examples and use them
for templates. Look up the syntax (Matlab
A common operation is to step through a number of cases and
remember the "best" or "worse" or "smallest", etc.
A standard approach is to create a "minimum" variable,
assign the first element in a list (vector, matrix...) to it, and
loop through the rest of the elements replacing the minimum if you
find a smaller element (Att 4.1.2).
function [min_val, min_index] = myminvec(vec)
min_val = vec(1);
min_index = 1; % Don't forget this!!
for index = 2:length(vec)
if vec(index) < min_val
min_val = vec(index);
min_index = index;
end % if
end % for
end % function
Again, matlab has
min
can take a matrix but then has
a possibly surprising result.
>> y = [ 4 -1
2 7];
>> min(y)
ans = [2 -1] % mins of columns (!)
>> min(min(y))
ans = -1
Add all positive elements of a vector vec
PosSum = 0;
for ndx = 1:length(vec)
if vec(ndx) > 0
PosSum = PosSum+vec(ndx);
end
end
Reverse a vector.
avec = [ 1 1 2 3 5 8 13];
N = length(avec);
backvec = zeros(1,N);
for forward_ndx = 1:N
backwards_ndx = N-forward_ndx+1;
backvec(backward_ndx) = avec(forward_ndx);
end
backvec
backvec =
13 8 5 3 2 1 1
The following example illustrates the perils of data-type conversion! What's supposed to happen here? What does happen here?? Very strange... (hint: Attaway 1.5.3.1: linear indexing)
for ndx_var = [1 2; 3 4]
ndx_var % just report its value
end
ndx_var =
1
3
ndx-var =
2
4
Nested FOR Loops (Attaway 4.2)
It is very common to nest for loops inside of
for loops e.g., to add all matrix elements
function sum = mat_add(A)
% Sum elements of A
% sum is name of MatLab builtin but no conflict (scoping!)
[NRows NCols] = size(A);
sum = 0; % initialize
for row = 1:NRows
for col = 1:NCols
sum = sum + A(row, col);
end % NCols
end % NRows
end % mat_add
This is same as the vectorized
sum(sum(A));
Note that nesting loops results in a dramatic increase in the number of operations carried out - the product of the the number of elements in the ranges. This gets out of hand very fast, so loops are seldom nested terribly deeply.
Please Don't Do This!!
Exercise 3.30:
function out = choose(in)
in = input('give me a number');
choice = menu('choose a function', 'ceil','round','sign');
.....
Functions talk to other functions! They don't need humans!
Diversity!
In Matlab:
for Mndxvar = 1:10
fprintf('\n Mndxvar = %d', Mndxvar);
Mndxvar = 2*Mndxvar;
end
yields
Mndxvar = 1
Mndxvar = 2
Mndxvar = 3
...
Mndxvar = 9
Mndxvar = 10
In C:
#include
yields
Cndxvar = 1
Cndxvar = 3
Cndxvar = 7
Why should we not be surprised?
The
The
Functions normally return when the code reaches the end of the function.
The
The
% script forit.m
for i = 1:5
disp(i*i);
end
i
fprintf('\n with break\n');
for i = 1:100
if (i*i) > 75
break
end
disp(i*i);
end
i % end of script
Running the script produces the following:
>> forit
1
4
9
16
25
i =
5
with break
1
4
9
16
25
36
49
64
i =
9
The
In the preceding example, the break destroys the implicit semantics of
the
The continue
statement is less
common than break
.
it aborts the current iteration of the loop and resumes iterating with
the next iteration of the loop.
Here we want to print out the indices of the positive values in a vector
vec = 2*rand(1,10) -1 % random numbers between -1 and 1
for i = 1:length(vec)
if vec(i)<0
continue
end
i % print if get here
end
Running the above produces:
vec =
Columns 1 through 6
0.52711 0.68345 0.48181 -0.47806
0.92415 -0.28688
Columns 7 through 10
-0.35907 0.62747 -0.64415 0.88629
i =
1
i =
2
i =
3
i =
5
i =
8
i =
10
Note that the code is a bit obscure and hard to follow, and this particular
program would be better written using an
This is often true, and if you find yourself wanting to use a
The most common "legitimate" use of these statements involves cleanly getting out of nested code when error conditions occur.
In contrast,
Since
The general form of a
while condition-expression
action
end
The condition must be true to get loop started and must become false sometime or you get infinite loop (exit one of these with CTRL-C if it happens to you in Matlab).
Specifically, the choice is:
Do N Times (where we know N in advance)
vs.
Do until job's done (and it's hard to say how long that will take)
Example: for vector X,
for i = 1:length(X)
X(i) = X(i) / 7;
end
divides every element of X by 7. Since we can easily know X's size,
a
On the other hand, suppose we want to add up 1/N for N =
1,2,3... etc. until the sum is only changing very slowly.
Stating this in terms of a change threshold we want to get below is
more natural than trying to figure out the appropriate N beforehand.
Hence a
denominator = 1;
term = 1/denominator;
sum = 0;
while term > .0001 % we're not done!
sum = sum + term;
denominator = denominator + 1;
term = 1/denominator;
end
sum
sum = 9.7875
(Aside: 1/1, 1/2, 1/3,... is the harmonic series from its place in the theory of music overtones (unison, octave, fifth, third...). What is its infinite sum?)
Other comments:
Work in
For this chapter and for the Pi project, please forget about "vectorizing" your code. Don't do it! Refer to matrix elements (Mat(i,j)), not matrices (Mat). And don't use the : range operator to refer to rows or columns of a matrix. Use loops. Don't create matrices for simple series like 'all odd numbers between 3 and 97' that can easily be done by repetition. On the other hand, for irregularly-spaced data like [ 1 2 5 10 20 50 100 200 1000], clearly a data vector is called for.
for x = 1:1000
for y = 1:1000
for z = 1:1000
statements; % A billion (x,y,z)s!!
end
end
end
We'll be seeing a lot of these, e.g. in Prog. Asst. 3.
Informally, statistics are numbers used to describe collections of other numbers, or data. Three common ones are the mean (average), the median, and the variance ( = standard_deviation^{2}).
We know about means and medians. Do we?. Variance and standard deviation measures how "spread out" the data are, how much the various data points differ from their mean value. If they're all the same, the variance is zero.
Naive implementation: read and use whole data set twice:
one 1:N
But! Notice, in the second moment equation above,
Var(X) = E[(X - μ)^{2}] = E[X^{2} - 2 μ X + μ^{2}] = E[X^{2}] - 2 μ E[X] + μ^{2} = E[X^{2}] - μ^{2} = E[X^{2}] - E[X]^{2}
Thus we can calculate mean and variance with only one pass through the data (one loop). In the loop that currently calculates the mean by summing X(i,j), we just need to accumulate another sum, this one sums the squares of the elements, X(i,j)*X(i,j). When done summing (after the for-loops), calculate the mean of the summed X's and X^{2}'s, do the subtraction and a square root, and we've got mean, variance, and stdev..
Some Theory for those interested
Almost always, one of these statements is best for what you want to do.
FLOWCHART FRAGMENTS
Destroy before reading or read and forget! This is very ugly stuff!
For some simple tests, Matlab has a vectorized version,
which unfortunately uses a form of indexing where matlab pretends
an n-dimensional matrix is really a long 1-D vector (matrix written
out columnwise). We've seen this weirdness before when we put a 2-d
matrix in as the 'range' in a
The following is for your horrified amusement only.
We do not recommend knowing about, let alone using, this
confusing nonsense.
>> a = magic(3)
a =
8 1 6
3 5 7
4 9 2 % magic square
>> g5 = a > 5 % employ matrix form of >
g5 =
1 0 1
0 0 1
0 1 0 % a binary matrix
% Treat relation as index!
>> a(a>5) % same as a(g5)
ans =
8
9
6
7 % returns values > 5 in a column
>> find(a>5) % returns indices of values >5
ans =
1
6
7
8
% BUT they're the 1-D form of indices!