#!/usr/bin/env python3
# -*- coding: utf-8 -*-
"""
Created on Tue Sep 10 14:10:38 2024

@author: sam
"""

def f(x):
    return x**3 + 1

def qsimpn(f,a,b,N):
    h = (b-a)/N
    I = f(a) + f(b)
    for i in range(1,int(N)):
        if i % 2 == 0: #for even points i
            xi = a + i*h #gets even point xi
            I += 2*f(xi) #adds 2 times the value of the function at xi
        elif i % 2 != 0: #for odd points i
            xi = a + i*h #gets odd point xi
            I += 4*f(xi) #adds 4 times the value of the function at xi
    I *= h/3 #multiplies I by h/3
    return I

result = qsimpn(f,0,2,1e6)
print(result)