Использование ODE в другом ODE с аргументами делает код очень медленным
Я написал следующий код. Это ODE, в котором есть параметр в качестве другого ODE.
Как мы видим, M(m0,z,b,c) используется в другом ODE, который сам является функцией ODE. код очень медленный, может ли кто-нибудь дать мне предложение, как его улучшить?
import numpy as np
from scipy.integrate import odeint
def model(m,z,c,b):
dmdt = ((c**2-m)/(1+z))*(6-9*(m/c**2)+3*b*(m+(m**2)))
return dmdt
def M(m0,z,c,b):
m = odeint(model,m0,[0,z], args= (c, b))
mm=m[-1,0]
return mm
def model1(H ,z,m0,c,b):
c = 0.6
b=0.035
dHdt = (H/(1+z))*(6-9*(M(m0,z,c,b)/c**2)+3*b*(M(m0,z,c,b)+(M(m0,z,c,b)**2)))
return dHdt
def model2(H0,z,m0,c,b):
H = odeint(model1,H0,[0,z], args=(m0,c,b))
HH=H[-1,0]
return HH
print(model2(70,1,0.75,0.69,0.035))
@Grisha Нет, в бесплатных параметрах больше 100000 номеров.
Если у вас есть один вызов дорогостоящей функции, вы не повторяете его с одними и теми же параметрами несколько раз, вы вызываете его один раз, сохраняете результат в переменной и повторно используете его там, где это необходимо. Но это все же не лучший способ ускорить вычисления, см. Ответ.
Кроме того, вы не должны изменять параметры в model1, параметры, переданные вызовом model2, становятся недействительными из-за жестко заданных значений.
Ответы 1
Вы можете решить связанную систему как связанную систему.
def model(U,z,c,b):
M, H = U
dMdt = ((c**2-M)/(1+z))*(6-9*(M/c**2)+3*b*(M+M**2))
dHdt = (H /(1+z))*(6-9*(M/c**2)+3*b*(M+M**2))
return [dMdt, dHdt]
def solution(H0,z,m0,c,b):
U = odeint(model,[m0,H0],[0,z], args=(c,b))[-1]
M, H = U
return H
print(solution(70,1,0.75,0.69,0.035))
который быстро возвращает 0.107569653042, в то время как ваш код с модификациями
def model1(H, z, m0, c, b):
mm = M(m0,z,c,b)
dHdt = (H/(1+z))*(6-9*(mm/c**2)+3*b*(mm+(mm)**2)))
return dHdt
возвращает аналогичный 0.107569746892 несколько медленнее. Эти 6 цифр совпадения соответствуют допускам ошибок по умолчанию 1e-6.
Для получения результатов с большей точностью установите контрольные параметры для допусков ошибок atol, rtol.
Для дальнейшего сокращения операций сделайте
def model(U,z,c,b):
M, H = U
factor = (6-9*M/c**2+3*b*(M+M**2))/(1+z)
return [(c**2-M)*factor, H*factor]
Если ваша задача действительно масштабная, используйте скомпилированный язык программирования для быстрого вычисления массовых чисел.
Спасибо. Я не мог уловить это предложение If your task is really massive, use a compiled programming language for fast mass number crunching.
Python хорош для быстрого прототипирования, в темные века для этой роли использовался BASIC. Но интерпретируемые языки, даже с JIT-компиляцией в байтовый код, имеют большие накладные расходы на каждую вычислительную операцию. Эти накладные расходы сводятся практически к нулю, когда вычисления компилируются в машинный код. Вы даже можете использовать тот же решатель, если можете работать с интерфейсом FORTRAN пакета решателей lsoda,.... Или используйте решатели gsl в C, обширной библиотеке решателей julia, ...
Ах, вы предлагаете использовать FORTRAN. Спасибо. Я делаю это для других расчетов.
Нет, вы можете использовать эти библиотеки с любым другим языком, вам нужно знать только о соглашениях о вызовах, особенно. выделение памяти или используйте оболочку, если она доступна, например scipy.integrate.odeint, которая является оболочкой python для lsoda. С другой стороны, в качестве теста я просто применил свой код к образцу 100x100 пар (b,c), и это заняло 24,6 секунды, поэтому, если вы выполняете массовые вычисления не очень часто, вы можете остаться с python и подождать 5-10 минут, которые потребуются. на 100000 проб.
Я могу понять, что вы говорите, но я не могу реализовать их своими собственными знаниями, я физик и просто использую то, что мне нужно, но из-за отсутствия надлежащих знаний я застрял в каком-то важном вопросе. Я не использую цикл for, я использую случайный выбор для каждого параметра и для 1000000 шагов. Я думаю, это огромно.
Это будет медленно - на каждом этапе решения большей модели вам придется решать меньшую. Сможете ли вы решить ее аналитически?