DZone Snippets is a public source code repository. Easily build up your personal collection of code snippets, categorize them with tags / keywords, and share them with the world

Snippets has posted 5883 posts at DZone. View Full User Profile

Simpson's Rule - Importing From CSV File

06.15.2008
| 3609 views |
  • submit to reddit
        
import csv

class simpsonCsv:
    def __init__(self, filename):
        reader = csv.reader(open(filename, "rb"))

        xlist = []
        ylist = []
        
        for row in reader:
            try:
               x = float(row[0])
               y = float(row[1])
            except TypeError:
                continue
            except ValueError:
                continue
            xlist.append(float(x))
            ylist.append(float(y))
        
        self.len = len(xlist) - 1
        self.xstart = xlist[0]
        self.xend = xlist[self.len]
        self.interval = self.xend / self.len
        
        self.data = map(None, xlist, ylist)
    
    def f(self, x):
        for (u,v) in self.data:
            if u == x:
                return v
                
        print "Bad x value (probably bad n): %s" % (x)
        return None
        
        
    
    def simpson(self, n):
        "Approximate the definite integral of f from a to b by Simpson's rule."

        if self.len % n != 0:
            print "Error: %d mod %d is not zero but should be." % (self.len, n)
            return -1
        
        
        h  = float(self.xend - self.xstart)/n
        
    
        si = 0.0
        sp = 0.0
        xk = 0.0
    
        for i in range(1, n, 2):
            xk = self.xstart + i*h
            si += self.f(xk)
    
        for i in range(2, n, 2):
            xk = self.xstart + i*h
            sp += self.f(xk)
        
        
        s = 2*sp + 4*si + self.f(self.xstart) + self.f(self.xend)

        return (h/3)*s
    

filename = "integral_tabela_CO2_csv.csv"
s = simpsonCsv(filename)

print s.simpson(20)
print s.simpson(40)
print s.simpson(60)
print s.simpson(100)
print s.simpson(300)
print s.simpson(600)
print s.simpson(900)
print s.simpson(1800)