cálculo de integrais em pyton

[Ignotum Peverel]

Como posso realizar um programa capaz de calcular uma integral pelo método de monte carlo?

Iniciei o programa dessa maneira:

 

from __future__ import print_function, division

import numpy as np

import matplotlib.pyplot as plt

 

%matplotlib inline

 

#imagem

a = -1

b = 1

 

#função

xr = np.linspace(a,b,1000)

yr = xr

 

c = np.min(yr)

d = np.max(xr)

 

n = int(input("Insira o número de pontos que deseja gerar aleatoriamente: "))

 

#número de rodadas

r = int(input("Insira o número de rodadas: "))

 

xa = np.random.random(n)*(b-a)+a

xb = np.random.random(n)*(d-c)+c

 

if #adicionar a condição para o cálculo da integral

 

print("O número total de pontos gerados é: ",(len(xa)))

print("O número de pontos contidos na área desejada é: ",(j))

 

plt.plot(xr,yr)

plt.plot(xa,xb,"+")

plt.show()

 

preciso fazer de modo com que o programa se repita um numero r de vezes para poder calcular a integral dada pela média das estimativas, o desvio padrão das estimativas, o desvio padrão da média para a estimativa, lembrando que a integral nada mais é que a área desejada da função. Como posso criar a condição "if" para que ele realize o cálculo e repita o programa r vezes?