2014 dxdy logo

Научный форум dxdy

Математика, Физика, Computer Science, Machine Learning, LaTeX, Механика и Техника, Химия,
Биология и Медицина, Экономика и Финансовая Математика, Гуманитарные науки




Начать новую тему Ответить на тему
 
 Проблемы с интерполяцией и поиском минимального значения
Сообщение24.11.2021, 15:25 


24/11/21
5
Всем доброго времени суток. Задача стоит в построении математической модели полимеризации полиизопрена и нахождения оптимального начального значения ДИБАГ (т.е s[:, 2]), который даст Среднемассовое значение (Mw) = 600000.
В теории начальное значение s[:, 2] (концентрация реагента) должна быть около 0,00378, чтобы получить Mw = 60000. На данном этапе, мой код, если его можно так назвать не показывает мне должного результата. По поиску минимального значения мне подсказали, вставил как в инструкции. Прошу помощи данного форума. Представляю сам код и вывод.
код: [ скачать ] [ спрятать ]
Используется синтаксис Python
import matplotlib.pyplot as plt
import numpy as np
from scipy.optimize import minimize
from scipy.integrate import odeint
from scipy.interpolate import splev, splrep, sproot
 
def  model_polisopren (s,t): # в данной функции записана Сист.Диф.Ур хим.процесса
    kp = 48; km = 0.0048; ka1 = 8.16; ka2 = 0.96 #константы скорости реакции
    dMdt = -s[0]*s[3]*(kp+km)-s[0]*(kp+km)*s[5]
    dA1dt = -ka1*s[1]*s[3]-ka1*s[1]*s[5]
    dA2dt = -ka2*s[2]*s[3]-ka2*s[2]*s[5]
    dPdt = -kp*s[0]*s[3]+(km*s[0]+ka1*s[1]+ka2*s[2])*s[5]
    dQdt = km*s[0]*s[3]+(ka1*s[1]+ka2*s[2])*s[3]
    dm0dt = kp*s[0]*s[3]-(km*s[0]+(ka1*s[1]+ka2*s[2]))*s[5]
    dn0dt = (km*s[0]+ka1*s[1]+ka2*s[2])*s[5]
    dm1dt = 2*kp*s[0]*s[3]+kp*s[0]*s[5]-(km*s[0]+(ka1*s[1]+ka2*s[2]))*s[7]
    dn1dt = (km*s[0]+(ka1*s[1]+ka2*s[2]))*s[7]
    dm2dt = 4*kp*s[0]*s[3]+kp*s[0]*s[5]+2*kp*s[0]*s[7]-(km*s[0]+(ka1*s[1]+ka2*s[2]))*s[9]
    dn2dt = (km*s[0]+(ka1*s[1]+ka2*s[2]))*s[9]
    return [dMdt, dA1dt, dA2dt, dPdt, dQdt, dm0dt, dn0dt, dm1dt, dn1dt, dm2dt, dn2dt]
   
s0 = [1.39,0.000177,0.00168,0.000007,0,0,0,0,0,0,0] # начальные значение концентрации реагентов
t = np.linspace(0, 5100, 10000)                     # задается время и шаг
s = odeint(model_polisopren, s0, t)                 # численное решение СДУ
Mw = ((s[:,9]+s[:,10])/(s[:,7]+s[:,8]+2**-100))*68  # считаем Mw тут
plt.plot(t,Mw,'b-',linewidth=2.0,label='Mw')
plt.xlabel("t")
plt.ylabel("Mw")
plt.legend()
plt.grid()
plt.show()
   
def poisk_min (s0):  # данная функция должна находить значения нач.концетрации при определенных значений
    x = splrep(t, Mw - 600000) # задаем условия, где нач.конц. должны быть такими, чтобы Mw = 60000 (мы считали его ранее)
    x0=sproot(x)
    y1=splrep(t, s[:, 2]) # указываем какая нач.концентрация нас интересует
    s1=splev(x0[0], y1)
    print (s1)
    return s1
res = minimize(poisk_min, s0, method='Nelder-Mead', tol=1e-3)
print (res)


И выводит следующие значения:

Код:
0.001676790666781049
final_simplex: (array([[1.39000000e+00, 1.77000000e-04, 1.68000000e-03, 7.00000000e-06,
        0.00000000e+00, 0.00000000e+00, 0.00000000e+00, 0.00000000e+00,
        0.00000000e+00, 0.00000000e+00, 0.00000000e+00],
       [1.39054297e+00, 1.77000000e-04, 1.68000000e-03, 7.00000000e-06,
        0.00000000e+00, 0.00000000e+00, 0.00000000e+00, 0.00000000e+00,
        0.00000000e+00, 0.00000000e+00, 0.00000000e+00],
       [1.39000000e+00, 1.77069141e-04, 1.68000000e-03, 7.00000000e-06,
        0.00000000e+00, 0.00000000e+00, 0.00000000e+00, 0.00000000e+00,
        0.00000000e+00, 0.00000000e+00, 0.00000000e+00],
       [1.39000000e+00, 1.77000000e-04, 1.68065625e-03, 7.00000000e-06,
        0.00000000e+00, 0.00000000e+00, 0.00000000e+00, 0.00000000e+00,
        0.00000000e+00, 0.00000000e+00, 0.00000000e+00],
       [1.39000000e+00, 1.77000000e-04, 1.68000000e-03, 7.00273437e-06,
        0.00000000e+00, 0.00000000e+00, 0.00000000e+00, 0.00000000e+00,
        0.00000000e+00, 0.00000000e+00, 0.00000000e+00],
       [1.39000000e+00, 1.77000000e-04, 1.68000000e-03, 7.00000000e-06,
        1.95312500e-06, 0.00000000e+00, 0.00000000e+00, 0.00000000e+00,
        0.00000000e+00, 0.00000000e+00, 0.00000000e+00],
       [1.39000000e+00, 1.77000000e-04, 1.68000000e-03, 7.00000000e-06,
        0.00000000e+00, 1.95312500e-06, 0.00000000e+00, 0.00000000e+00,
        0.00000000e+00, 0.00000000e+00, 0.00000000e+00],
       [1.39000000e+00, 1.77000000e-04, 1.68000000e-03, 7.00000000e-06,
        0.00000000e+00, 0.00000000e+00, 1.95312500e-06, 0.00000000e+00,
        0.00000000e+00, 0.00000000e+00, 0.00000000e+00],
       [1.39000000e+00, 1.77000000e-04, 1.68000000e-03, 7.00000000e-06,
        0.00000000e+00, 0.00000000e+00, 0.00000000e+00, 1.95312500e-06,
        0.00000000e+00, 0.00000000e+00, 0.00000000e+00],
       [1.39000000e+00, 1.77000000e-04, 1.68000000e-03, 7.00000000e-06,
        0.00000000e+00, 0.00000000e+00, 0.00000000e+00, 0.00000000e+00,
        1.95312500e-06, 0.00000000e+00, 0.00000000e+00],
       [1.39000000e+00, 1.77000000e-04, 1.68000000e-03, 7.00000000e-06,
        0.00000000e+00, 0.00000000e+00, 0.00000000e+00, 0.00000000e+00,
        0.00000000e+00, 1.95312500e-06, 0.00000000e+00],
       [1.39000000e+00, 1.77000000e-04, 1.68000000e-03, 7.00000000e-06,
        0.00000000e+00, 0.00000000e+00, 0.00000000e+00, 0.00000000e+00,
        0.00000000e+00, 0.00000000e+00, 1.95312500e-06]]), array([0.00167679, 0.00167679, 0.00167679, 0.00167679, 0.00167679,
       0.00167679, 0.00167679, 0.00167679, 0.00167679, 0.00167679,
       0.00167679, 0.00167679]))
           fun: 0.001676790666781049
       message: 'Optimization terminated successfully.'
          nfev: 103
           nit: 8
        status: 0
       success: True
             x: array([1.39e+00, 1.77e-04, 1.68e-03, 7.00e-06, 0.00e+00, 0.00e+00,
       0.00e+00, 0.00e+00, 0.00e+00, 0.00e+00, 0.00e+00])

 Профиль  
                  
 
 Re: Проблемы с интерполяцией и поиском минимального значения
Сообщение25.11.2021, 07:55 


21/05/16
4292
Аделаида
А вы уверены, что то, что ваша функция model_polisopren не использует t - правильно?

 Профиль  
                  
 
 Re: Проблемы с интерполяцией и поиском минимального значения
Сообщение25.11.2021, 09:03 


24/11/21
5
kotenok gav, да аргумент t используется правильно. В нем задается ось X с размерностью от 0 до 5100 с шагом 10000. (5100 это время, т.е секунды)

 Профиль  
                  
 
 Re: Проблемы с интерполяцией и поиском минимального значения
Сообщение26.11.2021, 08:11 


21/05/16
4292
Аделаида
Я не об этом. Он нигде не используется в теле функции.

 Профиль  
                  
 
 Re: Проблемы с интерполяцией и поиском минимального значения
Сообщение26.11.2021, 09:44 


24/11/21
5
kotenok gav, Аргумент t задается в Odeint для численного решения и затем для построения графика ))))

 Профиль  
                  
 
 Re: Проблемы с интерполяцией и поиском минимального значения
Сообщение26.11.2021, 13:14 
Заслуженный участник


04/03/09
906
Кажется, ваша минимизация делает совсем не то, что нужно. Уточните, требуется подобрать такое начальное значение s0[2] при неизменных остальных начальных значениях, чтобы Mw через 5100 секунд было как можно ближе к 600000?

 Профиль  
                  
 
 Re: Проблемы с интерполяцией и поиском минимального значения
Сообщение26.11.2021, 20:25 


24/11/21
5
12d3
Да, моя минимизация немного неправильная ((( А так, задача верная - нужно найти s0[2] при неизменных начальных значениях.

 Профиль  
                  
 
 Re: Проблемы с интерполяцией и поиском минимального значения
Сообщение27.11.2021, 00:02 
Заслуженный участник


04/03/09
906
Тогда, по сути, вам нужно численно решить уравнение $f(x) = 0$, где $f$ такова:
Используется синтаксис Python
def f(x):
    s0 = [1.39, 0.000177, 0.00168, 0.000007, 0, 0, 0, 0, 0, 0, 0]
    s0[2] = x
    t = np.linspace(0, 5100, 10000)              
    s = odeint(model_polisopren, s0, t)                
    Mw = ((s[:,9]+s[:,10])/(s[:,7]+s[:,8]+2**-100))*68
    return Mw[-1] - 600000
 

Для этого есть подходящая функция из пакета scypi: https://docs.scipy.org/doc/scipy/reference/generated/scipy.optimize.fsolve.html

 Профиль  
                  
 
 Re: Проблемы с интерполяцией и поиском минимального значения
Сообщение01.12.2021, 09:24 


24/11/21
5
12d3
Спасибо. Буду изучать )))) Вы не подскажите какой метод можно реализовать для двумерной минимизации ??

 Профиль  
                  
Показать сообщения за:  Поле сортировки  
Начать новую тему Ответить на тему  [ Сообщений: 9 ] 

Модераторы: Karan, Toucan, PAV, maxal, Супермодераторы



Кто сейчас на конференции

Сейчас этот форум просматривают: нет зарегистрированных пользователей


Вы не можете начинать темы
Вы не можете отвечать на сообщения
Вы не можете редактировать свои сообщения
Вы не можете удалять свои сообщения
Вы не можете добавлять вложения

Найти:
Powered by phpBB © 2000, 2002, 2005, 2007 phpBB Group