Расчет - букашка - PullRequest
Новая
       3

Расчет - букашка

3 голосов
/ 20 июля 2011

Я использую NumPy, чтобы сделать некоторые расчеты по нахождению перехвата Y, через Апертуру между большой коробкой и маленькой коробкой. У меня более 100 000 частиц в большой коробке и около 1000 в маленькой. И на это уходит много времени. Все self.YD, self.XD - это очень большие массивы, которые я умножаю.

PS: ind - это индексы значений, которые необходимо умножить. У меня было ненулевое условие перед этой строкой в ​​моем коде.

Есть идеи, как мне сделать этот расчет проще? 

YD_zero = self.oldYD[ind] - ((self.oldYD[ind]-self.YD[ind]) * self.oldXD[ind])/(self.oldXD[ind]-self.XD[ind])

Спасибо!

UPDATE

Будет использовать умножение, деление, вычитание и все такое Numpy. сделать это быстрее? Или, может быть, если я разделю расчет. например.

чтобы сделать это первым: 

    YD_zero = self.oldYD[ind] - ((self.oldYD[ind]-self.YD[ind])*self.oldXD[ind])

и тогда следующая строка будет: 

    YD_zero /= (self.oldXD[ind]-self.XD[ind])

Есть предложения?!

ОБНОВЛЕНИЕ 2

Я пытался выяснить это, сейчас, но не так много прогресса. Меня беспокоит то, что знаменатель: 

    self.oldXL[ind]-self.XL[ind] == 0

и я получаю странные результаты.

Другая вещь - ненулевая функция. Я проверял это некоторое время сейчас. Кто-нибудь может сказать мне, что это почти так же, как найти в Matlab

Ответы [ 3 ]

3 голосов
/ 20 июля 2011

Возможно, у меня неправильный конец палки, но в Numpy вы можете выполнять векторизованные вычисления. Удалите замкнутую петлю while и просто запустите эту ...

YD_zero = self.oldYD - ((self.oldYD - self.YD) * self.oldXD) / (self.oldXD - self.XD)

Это должно быть намного быстрее.

Обновление: Итеративный поиск корня по методу Ньютона-Рафсона ... 

unconverged_mask = np.abs(f(y_vals)) > CONVERGENCE_VALUE:
while np.any(unconverged_mask):
    y_vals[unconverged_mask] = y_vals[unconverged_mask] - f(y_vals[unconverged_mask]) / f_prime(y_vals[unconverged_mask])
    unconverged_mask = np.abs(f(y_vals)) > CONVERGENCE_VALUE:

Этот код является только иллюстративным, но он показывает, как вы можете применить итеративный процесс с использованием векторизованного кода к любой функции f, которую вы можете найти производной от f_prime. unconverged_mask означает, что результаты текущей итерации будут применяться только к тем значениям, которые еще не сходятся.

Обратите внимание, что в этом случае нет необходимости повторять, Ньютон-Рафсон даст вам правильный ответ на первой итерации, так как мы имеем дело с прямыми линиями. То, что у вас есть, является точным решением.

Второе обновление

Хорошо, значит, вы не используете Ньютона-Рафсона. Чтобы вычислить YD_zero (y-перехват) за один раз, вы можете использовать,

YD_zero = YD + (XD - X0) * dYdX

, где dYdX - градиент, который в вашем случае выглядит как

dYdX = (YD - oldYD) / (XD - oldXD)

Я предполагаю, XD и YD - текущие значения x, y частицы, oldXD и oldYD - предыдущие значения x, y частицы, а X0 - значение x диафрагма.

До сих пор не совсем понятно, почему вы должны выполнять итерацию по всем частицам, Numpy может выполнить расчет для всех частиц одновременно.

2 голосов
/ 20 июля 2011

Поскольку все вычисления выполняются поэлементно, должно быть легко переписать выражение в Cython.Это позволит избежать всех тех очень больших временных массивов, которые создаются при выполнении oldYD-YD и т. П.

Другая возможность - numexpr.

0 голосов
/ 21 июля 2011

Я бы определенно пошел на numexpr. Я не уверен, что numexpr может обрабатывать индексы, но держу пари, что (или что-то подобное) будет работать: 

import numexpr as ne

yold = self.oldYD[ind]
y = self.YD[ind]
xold = self.oldXD[ind]
x = self.XD[ind]
YD_zero = ne.evaluate("yold - ((yold - y) * xold)/(xold - x)")
...