El problema consiste en encontrar las raíces de una función usando el método de bisección. Este método no es el más eficiente para la tarea y justamente por ello lo usaremos, así se podrá apreciar aún más la mejora al momento de paralelizar las tareas.
Antes de continuar le recomiendo al lector que visite algunas páginas para que entienda el algoritmo, que es bastante simple.
Ahora bien, como función usaremos el coseno ya que como todos saben es una función periódica que corta el eje X (como dije antes, la idea es que "el problema" sea simple), como es costumbre lo desarrollaré en Python:
import math
def f(x):
return math.cos(x)
def biseccion(f, a, b, nIter, eps = 0.0001):
fa = f(a)
fb = f(b)
m = (a + b) / 2.0
fm = f(m)
n = 0
while n <= nIter and math.fabs(fm) > eps:
if fa * fm < 0:
b = m
fb = fm
else:
a = m
fa = fm
m = (a + b) / 2.0
fm = f(m)
n += 1
res = \
{
'r' : m,
'fr': fm,
'iteraciones': n,
'criterio': n > nIter and \
'iteraciones' or 'tolerancia'
}
return res
if __name__ == '__main__':
r = biseccion(f, 0, 4*math.pi, 100)
print r
Como podrán notar no está escrito de la manera más eficiente posible (hablando en términos de Python), sin embargo me interesa que sea fácil de entender. La función bisección recibe como parámetros la función a la cual le buscará las raíces, el intervalo de búsqueda [a, b], la cantidad máxima de iteraciones y por último la tolerancia, error o epsilon, como prefieran llamarlo, así tendremos ambos criterios de convergencia en caso de que el otro falle. La función retorna un diccionario con los valores obtenidos (la raíz, el valor de la raíz en el punto, el número de iteraciones y el criterio de convergencia.
Ahora bien, compliquemos un poco más el problema para hacerlo más retador para el cumputador. Supongamos que deseamos encontrar todas las raíces en un intervalo muy amplio, por ejemplo [-1000, 1000], el método de bisección solo nos va a encontrar una raíz ¿cómo buscamos las otras? Para responder a esto usamos el método de búsqueda incremental.
El método de búsqueda incremental, como su nombre lo indica, consiste en buscar haciendo pequeños incrementos por lo que al final nos quedará una serie de conjuntos más pequeños en los cuales hacer la búsqueda, esto no garantiza todas las raíces, sin embargo se aproxima bastante bien (al menos para el ejemplo que deseamos desarrollar).
El método nos quedaría como sigue:
def incremental(inicio, fin, nIntervalos, nIter):
inc = (fin - inicio) / nIntervalos
res = []
for i in range(nIntervalos):
r = biseccion(f, inicio, inicio + inc, nIter)
r.update({'a' : inicio, 'b' : inicio+inc})
res.append(r)
inicio += inc
return res
Recibe como parámetros el inicio y el fin del intervalo, el número de trozos en el que se va a cortar dicho intervalo y el número de iteraciones que se hará en cada búsqueda de las raíces. Nos retorna una lista de diccionarios con los resultados de cada ejecución del método de bisección, para las pruebas iniciales recomiendo usar parámetros más conservadores, pero hice una prueba de la siguiente manera:inicio = time.time()
r = incremental(-10000.0, 10000.0, 100000, 1000)
print "Tiempo transcurrido: ", time.time() - inicio, "s"
Aquí una captura de htop mientras se ejecuta la aplicación:
El computador en el que hice las pruebas es un Intel Core i3 con dos núcleos HT (por eso se ven 4) de los cuales únicamente se está usando uno a su máxima capacidad.
Se llevó aproximadamente 64.4822211266 s. Pues bien, ahora nos preguntamos: ¿hace falta ejecutarlo de manera secuencial? ¿no podría aprovechar de mejor manera varios núcleos?
La respuesta es un rotundo "sí", pero la respuesta al "¿cómo?" la veremos en el siguiente artículo.
