Tuesday, January 3, 2012

Fixed point iteration

A fixed point for a function is a point at which the value of the function does not change when the function is applied. More formally, x is a fixed point for a given function f if
and the fixed point iteration
converges to the a fixed point if f is continuous.
The following function implements the fixed point iteration algorithm:
from pylab import plot,show
from numpy import array,linspace,sqrt,sin
from numpy.linalg import norm

def fixedp(f,x0,tol=10e-5,maxiter=100):
 """ Fixed point algorithm """
 e = 1
 itr = 0
 xp = []
 while(e > tol and itr < maxiter):
  x = f(x0)      # fixed point equation
  e = norm(x0-x) # error at the current step
  x0 = x
  xp.append(x0)  # save the solution of the current step
  itr = itr + 1
 return x,xp
Let's find the fixed point of the square root funtion starting from x = 0.5 and plot the result
f = lambda x : sqrt(x)

x_start = .5
xf,xp = fixedp(f,x_start)

x = linspace(0,2,100)
y = f(x)
plot(x,y,xp,f(xp),'bo',
     x_start,f(x_start),'ro',xf,f(xf),'go',x,x,'k')
show()
The result of the program would appear as follows:
The red dot is the starting point, the blue ones are the sequence x_1,x_2,x_3,... and the green is the fixed point found.
In a similar way, we can compute the fixed point of function of multiple variables:
# 2 variables function
def g(x):
 x[0] = 1/4*(x[0]*x[0] + x[1]*x[1])
 x[1] = sin(x[0]+1)
 return array(x)

x,xf = fixedp(g,[0, 1])
print '   x =',x
print 'f(x) =',g(xf[len(xf)-1])
In this case g is a function of two variables and x is a vector, so the fixed point is a vector and the output is as follows:
   x = [ 0.          0.84147098]
f(x) = [ 0.          0.84147098]

Thursday, December 29, 2011

Book review: Numpy 1.5 Beginner's Guide

I got the chance to read the book NumPy 1.5 Beginner's Guide written by Ivan Idris and published by Packt Publishing last month. My impression of the book was quite positive. It's a book based on examples, which incrementally introduce all the main features of the library. It is written with a simple language, and it is always easy to understand.


Contents and structure

The organization and flow are good. The ten chapters contain well thought out examples that you can use as building blocks for your scientific computing projects. Every example is structured in this way:
  • An introduction to the problem that the example will solve.
  • The code, commented line by line.
  • The result of the code.
  • A short recap of how the problem has been solved.
  • And, sometimes, a multiple choice question to help the reader to test his own understanding.
There is no attempt at teaching the mathematics behind the examples. Every example is a "how to" that can help you to learn how to use the library and can save hours of searching through the official documentation and more complicated texts.


The chapters 1,2 and 3 contain the starting points to use NumPy. They explain how to install NumPy, how to handle the NumPy arrays and how to use some of the basic mathematical/statistical functions provided by the library. Chapters 4 through 7 cover the basics about handling matrices, how to load and write data, how to write universal functions and cover some of the basic modules that are discussed. Chapter 8 explains how to use the unit test functions provided by NumPy. Finally, chapters 9 and 10 (my favorites!) introduce how to integrate NumPy with Matplotlib and SciPy.

Who is this book for?

This book is aimed at people who know Python and need to start using scientific computing in their programs. It is also suitable for people who use another scientific computing environment, such as Matlab, and want quick-start introduction to NumPy.

Thursday, December 15, 2011

Polar Charts with matplolib

A polar system is a two-dimensional coordinate system, where there are two coordinates: the radial and the angular coordinates. The radial coordinate denotes the point distance from a central point (pole) and the angular coordinate denotes the angle required to reach the point from the 0 degree ray (polar axis). Let's see an example of how to make polar charts with matplotlib:
from pylab import figure,polar,show
from numpy import arange,pi,cos

theta = arange(0, 2, 1./180)*pi # angular coordinates
figure(1)
polar(3*theta, theta/5) # drawing a spiral
figure(2)
polar(theta, cos(4*theta)) # drawing the polar rose
show()
The result of this script consists of two charts. The first with the spiral
and the second with the polar rose