Несколько экземпляров scipy.integrate.ode

Question

Я хотел бы использовать экземпляры scipy.integrate.ode (или scipy.integrate.odeint) в нескольких потоках (по одному на каждое ядро ​​ЦП) для решения нескольких IVP одновременно. Однако в документации сказано: "Этот интегратор не является повторно входящим. Вы не можете иметь два экземпляра оды, использующих интегратор" vode "одновременно".

(Также odeint вызывает внутренние ошибки, если создается несколько раз, хотя в документации об этом не говорится.)

Есть идеи, что можно сделать?

Answer 1

Одним из вариантов является использование multiprocessing (т.е. используйте процессы вместо потоков). Вот пример, который использует map функция multiprocessing.Pool учебный класс.

Функция solve принимает набор начальных условий и возвращает решение, сгенерированное odeint, "Серийная" версия кода в основном разделе вызывает solve несколько раз, один раз для каждого набора начальных условий в ics, "Многопроцессорная" версия использует map функция multiprocessing.Pool экземпляр для запуска нескольких процессов одновременно, каждый вызов solve, map функция заботится о выдвижении аргументов solve,

Мой компьютер имеет четыре ядра, и по мере увеличения num_processes Максимальное ускорение составляет около 3,6.

from __future__ import division, print_function

import sys
import time
import multiprocessing as mp
import numpy as np
from scipy.integrate import odeint



def lorenz(q, t, sigma, rho, beta):
    x, y, z = q
    return [sigma*(y - x), x*(rho - z) - y, x*y - beta*z]


def solve(ic):
    t = np.linspace(0, 200, 801)
    sigma = 10.0
    rho = 28.0
    beta = 8/3
    sol = odeint(lorenz, ic, t, args=(sigma, rho, beta), rtol=1e-10, atol=1e-12)
    return sol


if __name__ == "__main__":
    ics = np.random.randn(100, 3)

    print("multiprocessing:", end='')
    tstart = time.time()
    num_processes = 5
    p = mp.Pool(num_processes)
    mp_solutions = p.map(solve, ics)
    tend = time.time()
    tmp = tend - tstart
    print(" %8.3f seconds" % tmp)

    print("serial:         ", end='')
    sys.stdout.flush()
    tstart = time.time()
    serial_solutions = [solve(ic) for ic in ics]
    tend = time.time()
    tserial = tend - tstart
    print(" %8.3f seconds" % tserial)

    print("num_processes = %i, speedup = %.2f" % (num_processes, tserial/tmp))

    check = [(sol1 == sol2).all()
             for sol1, sol2 in zip(serial_solutions, mp_solutions)]
    if not all(check):
        print("There was at least one discrepancy in the solutions.")

На моем компьютере вывод:

multiprocessing:    6.904 seconds
serial:            24.756 seconds
num_processes = 5, speedup = 3.59

Answer 2

SciPy.integrate.ode, кажется, использует решатели LLNL SUNDIALS, хотя SciPy не говорит об этом явно, но они должны, по моему мнению.

Текущая версия решателя одеодов CVODE 3.2.2 является реентерабельной, что означает, что она может использоваться для одновременного решения множества проблем. Соответствующая информация появляется в пользовательской документации для CVODE v3.2.0 (SUNDIALS v3.2.0).

Вся информация о состоянии, используемая cvode для решения данной проблемы, сохраняется в структуре, и указатель на эту структуру возвращается пользователю. В пакете cvode нет глобальных данных, и поэтому он является реентерабельным. Информация о состоянии, специфичная для линейного решателя, сохраняется в отдельной структуре, указатель на которую находится в структуре памяти cvode. Повторное появление cvode было вызвано ожидаемым расширением мультикомпьютера, но также важно в однопроцессорной среде, где две или более проблем решаются смешанными вызовами пакета из одной пользовательской программы.

Но я не знаю, поддерживает ли этот параллелизм SciPy.integrate.ode или другие решатели од, такие как scikits.odes.ode.