Posted: April 11th, 2008
def simpson(f, a, b, n):
"Approximate the definite integral of f from a to b by Simpson's rule."
if n % 2 != 0:
print "Ups: n must be even!"
return -1
h = (float(b) - a)/n
si = 0.0
sp = 0.0
for i in range(1, n, 2):
xk = a + i*h
si += f(xk)
for i in range(2, n, 2):
xk = a + i*h
sp += f(xk)
s = 2*sp + 4*si + f(a) + f(b)
return (h/3)*s
def f(x):
return x**4
ni = 50
nf = 1000000
n = ni
a = -20
b = 0
s = []
qc = []
ec = []
t = 0.0
div = 0.0
i = 0.0
while n < nf:
t = simpson(f, a, b, n)
s.append(t)
n = n * 2
for i in range(0, len(s)-2):
div = (s[i+2] - s[i+1])
if div == 0:
break
t = (s[i+1] - s[i]) / div
qc.append(t)
t = div / 15
ec.append(t)
for i in range(0, len(qc)):
print "%.12f, %.12f, %.12f => qc=%.12f, e=%.12f" % (s[i], s[i+1], s[i+2], qc[i], ec[i])
Posted: April 4th, 2008

My friend José Gaspar will be teaching an Advanced Administration Course (Curso Administração Avançada de Servidores Linux)
Some highlights/keywords:
- Servers: WWW, FTP, DNS, DHCP, LDAP, E-Mail
- Samba (with quotas)
- Security
- XEN
Location: Rua da Boavista, Porto Portugal (link)
Date: From April 22 to June 3.
Schedule: Tuesday and Thursday, 19h-23h (4h)
Duration: 72h
Price: 300 euros
More information (in portuguese):
http://moodle.libhertz.com
http://www.solutionsout.com/cursos.htm
Posted: April 4th, 2008
# this requires numpy get it from http://numpy.sf.net
from copy import deepcopy
from numpy import *
# this function, swapRows, was adapted from
# Numerical Methods Engineering with Python, Jean Kiusalaas
def swapRows(v,i,j):
"""Swaps rows i and j of vector or matrix [v]."""
if len(v) == 1:
v[i],v[j] = v[j],v[i]
else:
temp = v[i].copy()
v[i] = v[j]
v[j] = temp
def pivoting(a, b):
"""changes matrix A by pivoting"""
n = len(b)
for k in range(0, n-1):
p = int(argmax(abs(a[k:n, k]))) + k
if (p != k):
swapRows(b, k, p)
swapRows(a,k,p)
def gauss(a, b, t=1.0e-9, verbose=False):
""" Solves [a|b] by gauss elimination"""
n = len(b)
# make copies of a and b so as not to change the values in the arguments
tempa = deepcopy(a)
tempb = deepcopy(b)
# check if matrix is singular
if abs(linalg.det(tempa)) < t:
return -1
pivoting(tempa, tempb)
for k in range(0,n-1):
for i in range(k+1, n):
if tempa[i,k] != 0.0:
m = tempa[i,k]/tempa[k,k]
if verbose:
print "m =", m
tempa[i,k+1:n] = tempa[i,k+1:n] - m * tempa[k,k+1:n]
tempb[i] = tempb[i] - m * tempb[k]
# Back substitution
for k in range(n-1,-1,-1):
tempb[k] = (tempb[k] - dot(tempa[k,k+1:n], tempb[k+1:n]))/tempa[k,k]
return tempb
def residue(a, b, c):
"""Calculates the residue of a system solved by gauss elimination"""
n = len(b)
t = a * c # t is the A with the values of x replaced (an [n x n] matrix)
s = []
for i in range(0, n):
s.append(sum(t[i])) # s is the solution
res = b - s # res is the residue
return res
#a = array([[1.0, 2.0, 0.0],[-1.0, 2.0, 3.0],[1.0, 4.0, 1.0]])
#b = array([3.0, -1.0, 4.0])
#a = array([[-1.414214, 2, 0],[1, -1.414214, 1], [0, 2, -1.414214]])
#b = array([1.0,1.0,1.0])
#a = array([[2.0, 2.0, 2.0],[1.0, 1.0, 5.0], [2.0, 5.0, 1.0]])
#b = array([6.0, 7.0, 8.0])
a = array([[1.001, 2.001, 3.001],[0.999, 2.0, 2.999], [1.002, 1.999, 2.999]])
b = array([4.003, 4.001, 3.999])
x = gauss(a, b)
print "Solution = ", x
#sol = linalg.solve(a, b)
#print "linalg Solution = ", sol
y = residue(a, b, x)
print "Residue = ", y
u = gauss(a, y)
print "Residue destribution = ", u
z = gauss(a, b+y)
print "New Solution (with added residue) = ", z
y2 = residue(a, b+y, z)
print "Residule of new solution = ", y2
if linalg.norm(y2) < linalg.norm(y):
print "New solution has a smaller residue."
else:
print "Original solution has a smaller residue."
Posted: March 31st, 2008
Dave Winer recently broke the news that Google will soon introduce a “Google Web Services”, a competitor to Amazon Web Services. This isn’t much of a shock to anyone. The extra bit that is somewhat of a surprise is the price: free. He then explains that the reason it will be free is that they can pay for it in the reduced cost of integrating new acquisitions into their infrastructures which would become, effectively zero.
I don’t agree with that. At least I don’t think that is the main reason. I think it’s just the same old business model – basic service is free, premium is paid. The same model you see in Google Apps and a lot of other business on the web.
The basic service will be limited, in terms of traffic and/or storage (file and DB) etc. This will be the most widely used by startups. The premium service will just be the same pricing model amazon – pay as you go, the only different is that you don’t start paying at zero usage but rather at a higher threshold. When a startup becomes successful its needs will grow exponentially. They will be using more than the maximum allowed for free, and they’ll need to pay. But that’s alright because now that they are successful they can afford it. And the 1% or so that will be successful will have subsidized everyone that didn’t make it and still give Google a nice profit. Let’s not forget that Google’s infrastructure is already here and even tens of thousands of failed ventures won’t make a dent in it. Successful ventures will generate enough cash to upgrade that same infrastructure. Specially as the cost of hardware continues to drop.
Off course I’m sure people at Google also thought of the acquisition factor. Past acquisitions have taken a lot of time to be integrated. I think it took something like a year for writely to become “Google Docs“. That’s a very long time on the web.
Combine the free GWS with the powerful web development frameworks like Ruby on Rails that allow single individuals to create useful applications quickly and the new marketting oppurtunities that the web2.0 has created and the cost of trying won’t be measured in millions, thousands or even hundreds of dollars. It will be measured in terms of hours – the hours you “wasted” trying. And that, in many cases, won’t even be “time wasted” but rather “experience gained”.
It’s a brave new world indeed.
Posted: March 31st, 2008
JabberLogBot is a jabber bot that records messages sent to it in a database.
The idea came from Nuno Dantas who wanted a jabber bot to record quick notes. He talked about it at a Prt.Sc dinner last Wennesday.
There’s also a simple PHP file in there that displays the data. That file is currently just for show as it is rather “plain”.

UPDATED: added a screenshot of the web viewer.
Posted: March 28th, 2008
The second beta version of the iPhone SDK is now available and includes Interface Builder, a powerful tool that allows you to visually build your interface and makes creating a UI as simple as drag and drop.
Screenshots:

Download at iPhone Dev Center
Posted: March 28th, 2008

Posted: March 27th, 2008

I have 6 evernote (previously mentioned here) invites up for grabs. If you’re interested, drop a comment with your email.
Posted: March 27th, 2008

21h00, Saturday, March 29
Duration: 2h30min
Prices:
Plateia – 30,00€
Tribuna – 25,00€
Camarotes 1a – 25,00€
Frisas – 22,50€
Galeria/Geral – 15,0
More Info @ Coliseu do Porto
www.kickboxingcup.com.pt
Posted: March 27th, 2008
I finally decided to take the Numerical Methods (MNUM) course. It turns out it’s a lot more fun than I thought. There is programming involved but you can chose to use whatever language you want. This is yet another nice excuse for me to use Python instead of C++ or Java. Last semester I was able to use Python to implement the game logic for Software Application Laboratory (LAS), which is mostly an OpenGL course with IPC via sockets thrown into the mix, and to write an article on dynamic languages (focusing mostly on Python) for Software Engineering (ESOF).
But back to this semester, 3 classes into the semester and the teacher is already said something like “I’m going to learn python now. I didn’t believe when I heard someone saying it was the best language in the world, but now I see there might be some truth to that claim”. That and I suspect his next laptop might be a macbook but that’s another story.
There are a few things that make Python great for Numerical Methods. In my opinion, Python’s clear, easy to understand, syntax is the most important one.It makes algorithms easier to implement. The syntax ends up being very close to language neutral pseudocode available in numerical methods books. Also Python’s datatypes as well as those provided by other libraries can be very useful.
The following code implements the stuff in chapter 2 (determining zeros) of the course. The methods implemented are Bisection, Rope and Newton. The function returns both the solution and the number of iterations necessary to get to that solution.
UPDATE: forgot the book - Numerical Methods in Engineering with Python
Appendix A – mnum2.py
from math import log
def bisect(f, a, b, e):
""" Determines zero between a and b using Bisection. """
n = 0
fa = f(a)
if fa == 0.0: return (a, n)
fb = f(b)
if fb == 0.0: return (b, n)
while (abs(a-b) > e):
c = 0.5*(a+b)
fc = f(c)
if fc == 0.0: return (c, n)
n = n + 1
if fb*fc < 0.0:
a = c
fa = fc
else:
b = c
fb = fc
if fa < fb:
return (a, n)
else:
return (b, n)
def rope(f, a, b, e):
""" Determines zero between a and b using the Rope methode. """
n = 0
fa = f(a)
if fa == 0.0: return (a, n)
fb = f(b)
if fb == 0.0: return (b, n)
while (abs(a-b) > e):
c = (a*fb - b*fa) / (fb - fa)
fc = f(c)
if fc == 0.0: return (c, n)
n = n + 1
if fb*fc < 0:
a = c
fa = fc
else:
b = c
fb = fc
if fa < fb:
return (a, n)
else:
return (b, n)
# Note: must verify that for the function f and guess c
# the method will _converge_.
def newton(f, df, c, t):
""" Determines zero between a and b using Newton """
n = 0
fc = f(c)
if fc == 0.0: return (c, n)
while (True):
fc = f(c)
dfc = df(c)
if dfc == 0:
print "dfc is 0"
return (0, -1)
dc = -fc/dfc
c = c + dc
n = n + 1
if abs(dc) < t: return (c, n)
##Tests
#def f(x): return -log(x)+4.0
#def df(x): return -1.0/x
#x= bisect(f, 1, 70, 0.00000001)
#print x
#x = rope(f, 1, 70, 0.00000001)
#print x
#x = newton(f, df, 0.1, 0.0001)
#print x
Posted: March 26th, 2008
Ok, from the comments on the previous post it seems people are NOT “getting it”. This is NOT a “problem”, this is funny.
So I guess people don’t know what the functions min() and max() do or are just confusing them
max([...]) - this function usually takes a list of numbers and returns the highest number in the list.
e.g.
if x belongs to [0,20] then y = max(x, 18) => y belongs to [18-20]
max(20,18) = 20
max(19, 18) = 19
max(10, 18) = 18
min([...]) - this function usually takes a list of numbers and returns the lowest number in the list.
e.g.
if x belongs to [0,20] then y = min(x, 18) => y belongs to [0-18]
min(20,18) = 18 [<- fixed thanks Mind Booster Noori]
min(19, 18) = 18
min(10, 18) = 10
This means that according to the function the function displayed in my previous post, according to that function, ALL students grades will be between 18 and 20. That means if you have a freaking 0 in the exam and a 0 in the Assignment you’ll get an 18/20 (which is an excellent grade). i.e. THE MINIMUM grade is 18. It’s a simple mistake in the equation. What the page should say to be correct is
Nota Final = min(Nota do exame final+Nota do trabalho, 18)<- this would be the correct equation
And hell this was posted under “Entertainment”.
So the previous post wasn’t supposed to be about something “unfair” it was supposed to be about making fun of a simple mistake that under the circumstances really is funny.
Posted: March 25th, 2008
Or not. Someone got min() and max() confused:

Nota Final = max(Nota do exame final+Nota do trabalho, 18)
translation: Final Grade = max(Exam Grade + Assignment Grade, 18)
[Note: it’s 18/20 – 20 is the maximum grade possible)
Than again, maybe it’s not a mistake. Perhaps the teacher really enjoys giving high grades in the 18-20 range. Now wouldn’t that be different…
Link to the page
UPDATE: Read my next post for an explanation of why this is funny.
Posted: March 25th, 2008
I just got an email (via the faculty-wide spam network) about the ACM International Collegiate Programming Contest (ICP) or more precisely about SWERC 08. For a few moments I thought “hey, this might be fun”. A second later I read something like
“Languages allowed: C, C++, Java or Pascal”
Yay it’s 1995 again!!!
Pascal? Didn’t Pascal die like more than a decade ago?
C++ and Java are so 90′s… where are the cool languages? Python? Ruby? Heck even Erlang and Haskell?
Considering the “retro” theme, I would’ve was expecting Lisp and Fortran to be in there. (Not (that I have anything against (lisp))).
So basically the point of the contest is seeing how well you do with a handicap? Kinda like the Olympic 100m with a broken leg (or 2 if you’re using C)?
Sounds very interesting but I rather go cut myself with the portuguese tokio hotel fans (via Tiago Farrajota)…
Posted: March 24th, 2008
If you happen to study at FEUP or maybe even somewhere else in Portugal or you just happen to enjoy this kind of humor check out:
Motivational Posters (warning1: not work safe) (warning2: mostly in portuguese)


