Lecture 5: Procedures for performing iterative calculations.

Procedural languages like Fortran, C and C++ use Do, For, While. Mathematica provides these function, but the syntax may be a bit different.

The Do loop

Do[expression,{i,1,n}]
evaluates expression for i=1, then for i=2, .... then for i=n.
It is easier to understand this by looking at examples.

Do[Print[i],{i,1,6}]

Do is often used to calculate the sum of the elements of a list

li={2,4,6,8,10};
s=0;
Do[s=s+li[[i]],{i, 1 , Length[li]}]

Here we started by initializing s (which will contain the sum, in the end) to s = 0. When i = 1 the computer evaluates
	s = s +li[[1]] which gives s =0 + 2 = 2;
Next the computer makes i = 2 and evaluates
	s = s + li[[2]] = 2 + 4 = 6
Next it makes i =3 and evaluates
	s = s + li[[3]] = 6 + 6 = 12
etc.  Final result is that s = 0 + 2 + 4 + 6 + 8 + 10 = 30. Indeed, this is the case as you can see below.

s

Note that the Do loop has no output. We assign values to s (inside the Do loop) and then we use s when we need it.

It is easier, faster, more readable and more elegant to use what is called functional programing, which performs the sum as shown below:

Plus@@li

A Do loop is often use to multiply the elements of two lists and add the result of the multiplications:  I start with the lists

 	li={2,4,6,8,10}
 	
 and
 
 	li2={a,b,c,d,e}
 	
 and I want to calculate 
 
 	w= 2 a+4 b+6 c+8 d+10 e
 	
 Here is how you do this with a Do loop:

li={2,4,6,8,10};
li2={a,b,c,d,e};
s=0;
Do[s=s+li[[i]]*li2[[i]],{i,1,Length[li]}]

s

The functional approach is again simpler, more elegant and faster:

li.li2

Another use of a Do loop is to start with two lists

	li={2,4,6,8,10}
 and
 	li2={a,b,c,d,e}
 	
 and  create a new list containing
 
 	w={ 2 a, 4 b, 6 c, 8 d, 10 e}

First I create a list w that has the proper length and give arbitrary values to its elements (here I make all elements in w equal to zero).  Then I give the elements of w the values 2 a, 4 b, etc. 
I'll tell you soon why this is a good strategy when dealing with large lists.

w=Table[0,{i,1,Length[li]}]
Do[w[[i]]=li[[i]]*li2[[i]],{i,1,Length[li]}]
w

There is no output from the Do loop but we have captured what we want in w.

Again, this is easier done by functional programing:

li*li2

Another way of doing the above calculation is to use:

Append[{1,2,3},4]

p={};
li
li2
Do[p=Append[p,li[[i]]*li2[[i]]],{i,1,Length[li]}]
p

The trouble with this is that Append[list,a] copies the list and then adds a to it. This coping takes time when the list is long. The procedure we used first is more efficient.

You can have multiple Do loops and they are used, for example, to multiply a matrix  with a vector. For example: the matrix is 

	li={{x,y,z},{a,b,c},{A,B,C}}
	
and the vector is

	vec={M,N,K}
	
and the product 

	li . vec={M x + N y + K z, c K + a M + b N,  C K + A M + B N} 
	
Here is how you do this in a Do loop

li={{x,y,z},{a,b,c},{A,B,C}};
li//MatrixForm
vec={M,N,K};
vec//MatrixForm
w={0,0,0};
Do[w[[i]]= w[[i]]+li[[i,j]]*vec[[j]],{i,1,3},{j,1,3}] 
w//MatrixForm

It is a lot simpler to do this in a functional programing.

li.vec//MatrixForm

The For loop

Another procedure which is obsolete in Mathematica but is alive and well in C and C++ is the For loop. The syntax is:

	For[start, test, increment, expression]
	
An example is

s=0;
For[i=3,i<6,i=i+1,Print["i = ",i," first s= ",s];s=s+i^2; Print["i= ",i," s= ",s]]

s

To see what happened I place some Print statements strategically inside the loop. The loop starts by making i = 3 as we instructed it in For[i = 3, test, increment, expression]. It tests and finds that i is less than 6 and accepts this value of i. Then it goes and evaluates s = s+i^2  = 0 + 3^2 = 9;  Next it goes and increments i to i + 1 (which is now 4). It tests that 4 < 6 and accepts i = 4. Then it goes and evaluates s=s+(i^(2  ))= 9+(4^(2 ))= 9 + 16 = 25, then it increments i to i+1 etc. The loop stops when i = 6 because the test 6<6 fails.

You can see this better in the following loop

For[k=2,k<=5,k=k+1,Print["k= ",k]]

Remember that k only takes values from the initial one (k=2 in our example) to the last one that gives True when it is tested.

The increment is up to you. You can use

For[k=-2,k<20,k=k+3,Print["k= ",k]]

Like in the case of Do you could nest For loops

For[k=1,k<5,k=k+1,For[i=2,i<7,i=i+2,Print["{k,i}= ","{",k," , ",i," }"]]]

Here is how you multiply a matrix with a vector by using For
 For example: the matrix is 
 
	li={{x,y,z},{a,b,c},{A,B,C}}
	
and the vector is

	vec={M,N,K}
	
and the product 

	li . vec={M x+N y+K z,c K+a M+b N,C K+A M+B N} 
	
Here is how you do this in a "double" For loop

li={{x,y,z},{a,b,c},{A,B,C}};
vec={M,N,K};
w={0,0,0};
For[i=1,i≤3, i=i+1,For[j=1,j≤3,j=j+1,w[[i]]= w[[i]]+li[[i,j]]*vec[[j]]   ]]
w

Here is an example of using the For loop in PChem.  You are working with the reaction 

	(Subscript[C, 2])(Subscript[H, 6])+(7/2)(Subscript[O, 2 ])→2 (Subscript[CO, 2]) +3(Subscript[H, 2])O
	
You want to write a function that given the list of stoichiometric coefficients will print out the reaction as seen above.  You want this function to be general and work for any reaction.

comp={"C2H6","O2","CO2","H2O"}

ν={-1,-7/2,2,3}

The function is given below. I will give you the basic idea of the program. I can use the sign of the stoichiometric coefficients to find whether a compound is reactant or product. I put the negative stoichiometric coefficients in a list Lν, the corresponding reactants in the list Rν; the positive stoichiometric coefficients and the corresponding products are placed in the lists Lν and Rν. After this sorting,  I print -Lν.Lcomp and Rν.Rcomp and put an arrow between them. That prints the text that I want. The For loop is used to let me run through the lists ν and comp and sort them.

reaction[ν_List,comp_List]:=Module[{Lν,Rν,Lcomp,Rcomp},
Lν={};Rν={};Lcomp={};Rcomp={};
For[i=1,i≤Length[ν],i=i+1,If[ν[[i]]<0,Lν=Append[Lν,ν[[i]]];Lcomp=Append[Lcomp,comp[[i]]],Rν=Append[Rν,ν[[i]]];Rcomp=Append[Rcomp,comp[[i]]]]] ;
Print[-Lν.Lcomp," → ", Rν.Rcomp]]

reaction[ν,comp]

Exercise: write a function doing the same thing as reaction, but use Select to do the sorting.

While

Another function that allows iterative calculations is While. 
This has the syntax:

	While[test, expression]
	
The expression is evaluated for as long as the test returns False.  
While  is used very much like the For loop, except that we initialize the loop outside of the While function (see i=4 in the example below), and one increments the iterator in the expression. The example below calculates 

	(a^4)+(a^7)+(a^10)+(a^13)+(a^16)+(a^19)

i=4;s=0;
While[i<20, s=s+a^i ;i=i+3; Print["s= ",s , " i= ",i]]

s

You can do with while pretty much what one can do with For, including nesting them. For example, here is how you calculate
	"s= "a (b^3)+(a^3) (b^3)+a (b^6)+(a^3) (b^6)+a (b^9)+(a^3) (b^9)

i=1; j=3;Clear[a,b,s];s=0;While[i<5,
 j=3;While[j<10,s=s+a^i*b^j;j=j+3;Print["i= ",i,"; j= ",j,"; s= ",s]];i=i+2;]

Nest and NestList

Often in Mathematics and Numerical Analysis we need to evaluate expressions of the form 

	f[f[f[f[x]]]]
	
where f is applied successively n times (in the example above n = 4).
You may thing that this is a rarely used operation, but as you learn more mathematics and physics you will often use it. Mathematica provides a function, called Nest, for performing it. Moreover, if you want the list

	{x, f[x], f[f[x]], f[f[f[x]]], f[f[f[f[x]]]], ...}
	
you can get it by using the function NestList. Here is how they work

Nest[f,x,4]

%/.f->Sin

NestList[f,x,4]

Here are a few simple examples showing how these functions work:

Nest[#^2&,x,3]

NestList[#^2&,x,3]

Here is an application: the   sequence

	(Subscript[x, n+1]) =(1/2)((Subscript[x, n])+(b/(Subscript[x, n])))
	
started with 

	(Subscript[x, 0 ])= 1
	
converges to (Sqrt[b]). This provides an efficient method for calculating the square root on a computer. The function below has two arguments. b is the number whose square root we want to calculate and n is the number of iterations we want to use.

sqrt[b_,n_,x0_]:=NestList[(#+b/#)/2&,x0, n]//N

sqrt[4,5,1.]

sqrt[10,8,3]

There are other useful functions performing iterations such as  NestWhile, Fold, FoldList and FixedPoint. You can look them up in Mathematica book.