In [1]:
import numpy as np
import matplotlib.pyplot as plt
%pylab inline


Populating the interactive namespace from numpy and matplotlib

In [18]:
def funcion(a,b,c,d,e,f,g,h,i,j):
    return (a+b+c+d+e+f+g+h+i+j)**3

In [58]:
def Integra_10D(fun, lim1, lim2, NPoints, intentos):
    resultados = []
    for i in xrange(intentos):
        a = np.random.uniform(lim1,lim2,NPoints)
        b = np.random.uniform(lim1,lim2,NPoints)
        c = np.random.uniform(lim1,lim2,NPoints)
        d = np.random.uniform(lim1,lim2,NPoints)
        e = np.random.uniform(lim1,lim2,NPoints)
        f = np.random.uniform(lim1,lim2,NPoints)
        g = np.random.uniform(lim1,lim2,NPoints)
        h = np.random.uniform(lim1,lim2,NPoints)
        i = np.random.uniform(lim1,lim2,NPoints)
        j = np.random.uniform(lim1,lim2,NPoints)
        resultados.append((lim2-lim1)**10/NPoints*sum(fun(a,b,c,d,e,f,g,h,i,j)))
    return mean(resultados)

In [134]:
N = 10000
intentos = 20
res = Integra_10D(funcion,0.0,2.0,N,intentos)
print ('Valor de la integral con {} puntos y {} intentos = {}'.format(N,intentos,res))


Valor de la integral con 10000 puntos y 20 intentos = 1124465.76448

In [122]:
eNes= []
valores = []
for i in xrange(1,14):
    n = 2**i
    valores.append(Integra_10D(funcion, 0.0, 2.0, n, intentos))
    eNes.append(n)

In [128]:
plt.plot(eNes, valores)
plt.title(u'Valor integral vs. Número de puntos')
plt.xlabel('N')
plt.ylabel('Valor I')
plt.savefig('num_integral.pdf')
plt.close()



In [131]:
ValorAnalitico = 1126400.0
errores = 100*abs(ValorAnalitico - np.array(valores))/ValorAnalitico
x = 1/sqrt(np.array(eNes))

In [135]:
plt.plot(x, errores)
plt.title('Errores porcentuales vs. $1/\sqrt{N}$')
plt.xlabel('$1/\sqrt{N}$')
plt.ylabel('Error %')
plt.savefig('err_integral.pdf')
plt.close()

In [ ]:


In [ ]: