Отслеживайте успешные шаги для функции resolve_ivp в scipy.integrate - PullRequest
0 голосов
/ 25 мая 2020

У меня проблема, требующая решения системы ODE, и я использую для нее scipy.integrate.solve_ivp. При вызове интегратора я указываю интервал интегрирования, а также моменты времени, в которые должно быть оценено решение. Вызов выглядит примерно так:

Z_solution  = solve_ivp( fun= lambda t, y:  rhs_of_differential_equation( t, y, args_tuple  ),
                         t_span= time_interval, y0 = initial_array, 
                         method='RK45',
                         t_eval = np.array( [ 0.0 ] ),                                               
                         max_step   = parameters.max_step, 
                         first_step = parameters.initial_step,                                                                  
                         rtol = parameters.adaptive_step_relative_tolerance, 
                         atol = parameters.adaptive_step_absolute_tolerance,
                         )

Это означает, что меня интересует только решение при t = 0, а все остальное между ними не имеет отношения к моему случаю. Тем не менее, при вызове метода resolve_ivp я отслеживаю, какие значения t использует интегратор, поэтому я отслеживаю все «промежуточные» шаги в алгоритме. Но я хотел бы знать, есть ли способ действительно увидеть, какие шаги были успешными, а не промежуточные. Мне кажется, что приведенный выше вызов интегрирует все до t = 0.0, но не говорит вам, был ли этап интеграции успешным или нет ... Я был бы признателен за некоторую помощь здесь, спасибо!

ДОПОЛНИТЕЛЬНАЯ ИНФОРМАЦИЯ : Мой вектор решений y фактически зависит от сетки, которая также масштабируется с t. То есть, если y - мое решение ODE, оно будет состоять из N точек, где N - размер сетки (для фиксированного значения t). Это означает, что для конечной точки интегрирования t = 0,0 мое решение y представляет собой массив из N элементов. Но сетка, в которой вычисляются значения, также зависит от значения t. Есть ли способ последовательно изменять сетку на протяжении всего алгоритма интеграции? Я попытался передать сетку в массиве начальных условий, но, похоже, не работает должным образом ..

EDIT-2: Здесь я написал небольшой код, чтобы понять это. Моя проблема имеет очень большие изменяющиеся масштабы, и поэтому я хочу позволить интегратору выбирать необходимые временные шаги для меня, поэтому я НЕ хочу заранее выбирать t_evals для фиксированного числа точек; Мне нужно только начальное и конечное время интервала. В промежутках между ними интегратор будет выполнять необходимую ступенчатую регулировку, чтобы не превышать допуски. Затем моя проблема хочет идентифицировать эти внутренние точки (точки между [0,10], которые интегратор смог вычислить). Кажется, что объект SolveODE внутри связки, возвращаемой функцией solution_ivp, имеет значения "ts", а также метод call , который в принципе извлекает сгенерированное решение на этих временных шагах. Вот код для простой функции экспоненциального затухания:

import numpy as np
from scipy.integrate import solve_ivp
import matplotlib.pyplot as pyt


def my_fun2( t, y, *args ):    
    return -0.5*y    

vals = np.array( [ 100 ] )    
args_tuple = ( 0.01, 0.5 )

Z_solution  = solve_ivp( fun= lambda t, y: my_fun2( t, y, args_tuple  ),
                                         t_span= [0.0, 10], y0 = vals, 
                                         method='RK45',
                                         t_eval =   [0.0, 10],
                                         dense_output  = True,                                                                  
                                         rtol=1e-7, atol = 1e-7,                                           
                   )  

OdeSol  = Z_solution.sol
y_sol = OdeSol.__call__( OdeSol.ts  )[0] 
ts    = OdeSol.ts 

x = np.linspace( 0, 10, 1000 )
analytic = 100*np.exp( -0.5*x )

for s in range( 0, len(vals) ):    
    pyt.plot( ts, y_sol, 'ob', x, analytic, '-r' )
    pyt.show()

и решение аналити c и точки по (ts, y_sol) дают:

график кода

Тогда, если я правильно понял, синие точки - это точки, где интегратор нашел решения в интервале [0,10] БЕЗ превышения допусков. Если вы уменьшите допуски, то, конечно, количество синих точек будет уменьшено, может ли кто-нибудь подтвердить это?

Вдобавок ко всему, что мне нужно для моей системы, так это каким-то образом сообщать функции в каждой точке ts для остановки выполните некоторую корректировку сетки y_sol, найденного в этой точке, и продолжите интегрирование. Есть ли способ сделать это с помощью resolve_ivp? Я могу только думать о том, чтобы включить эту функциональность сам в resolve_ivp, но я действительно не хочу этого, так как это означало бы изменение кода resolve_ivp, чтобы обеспечить такую ​​функциональность, большое спасибо!

Добро пожаловать на сайт PullRequest, где вы можете задавать вопросы и получать ответы от других членов сообщества.
...