Explanation of New Concepts

Here we explain some new features in the code nr.py.

We use a for loop here to impose a maximum number of Newton Raphson iterations, and prevent the code from running on endlessly.




  for i in range(N):
    ...

This is a for loop. It specifies that the statement block indented just below the for statement (represented here by \ldots{}) should be executed repeatedly with i = 0, 1, 2, ..., N-1.

More generally, the for statement specifies a list of values that the variable (in this case i) should take. Here we use Python's range function to generate the list. With only one argument, range takes on integer values starting from 0 and ending one value short of the given end point. With two arguments, the first specifies the starting value, and with three arguments, the third specifies the step or increment. So, for example,




  range(5) = 0,1,2,3,4
  range(1,8) = 1,2,3,4,5,6,7
  range(1,9,2) = 1,3,5,7

We could accomplish the same thing with a while loop. Here is the logically equivalent pattern:

  i = 0    # initialization
  while i < N:
    ...
    i = i + 1   # incrementing
Notice that here we have to take care to place initialization and incrementing statements in the correct place. And it takes three statements to do the job of one for statement. For those reasons the for loop is preferable here.

The generic form is

  while(logical-expression)
    statement-block
The while loop repeats the statements inside the statement block as long as the logical expression is true. It is up to the programmer to assure that there is an escape from the loop, so it doesn't run on indefinitely - another reason to be cautious about using a while loop where a for loop will do.




    p = pnew
    ...
    pnew = p - f/dfdx

Notice that the values of p and pnew are being replaced repeatedly with this algorithm. The idea is that pnew always contains the most up-do-date guess for the solution. We keep p, the previous guess, so we can calculated its difference from [new.




    if math.fabs(p-pnew) < tol:
      ...
    }

This test is called the stopping condition for the method. We stop when the absolute value of the change in the estimated root is less than the specified tolerance. We want the absolute value of the change, since it could be positive or negative. The math library has the absolute value fuction math.fabs.




      sys.exit(0)

This statement forces an exit from the main program. The value 0 is the exit code. It is common practice to use an exit code of 0 for a normal completion and a nonzero code for an abnormal condition. Notice that at the bottom of the code we give an exit code of 1 to signal that the method did not converge. You can check the exit code of your program in the shell tcsh using the Unix command
   echo $status
immediately after your program exits. In the bash or sh shell the command is echo $?.

The break statement provides another way out of a loop. We could have written

  for i in range(N)
    ...  
    if math.fabs(p-pnew) < tol):
       break
The break statement immediately quits the statement block controlled by a for and while, and proceeds with the next statement after the indented statement block.

The logical problem in this case is that we reach the same point in the code if the process doesn't converge. So we need to check the value of i to be sure that didn't happen. For example, we could end it this way:

  for i in range(N):
    ...  
    if math.fabs(p-pnew) < tol):
       break
    ...
  if i >= N:
      print(Failed to converge after" << N << "iterations.")
      sys.exit(1)
  print("Root is", pnew, "to within", tol)
  sys.exit(0)

For completeness we mention the continue statement. It provides another way to skip around in a loop. While break skips the rest of the statement block and quits the loop, the continue statement skips the rest of the statement block and proceeds with the next iteration. Here is an example of its use in the same code:

  for i in range(N):
    ...  
    if math.fabs(p-pnew) > tol:
      continue
    print("Root is", pnew, "to within", tol)
    sys.exit(0)
  print("Failed to converge after", N, "iterations.")
  sys.exit(1)
Notice that the inequality is reversed so the loop continues as long as the shift in the root is larger than the tolerance and i < N. I find the original coding style clearer, however, since the code handling the special case is set off clearly with indentation. The continue statement is more appropriate when the standard nesting of statement blocks is awkward or when continuation is the less likely alternative, which is not the case here.