Наименьший_квадрат неточен в химическом виде ⇐ Python

[Anonymous]

У меня возникла проблема с использованием функции наименьших_квадратов для решения системы нелинейных уравнений с целью получения химического состава (в данном случае системы никель-аммиак). Система уравнений приведена здесь:

Код: Выделить всё

def equations(x):
Ni0, Ni1, Ni2, Ni3, Ni4, Ni5, Ni6, NH3 = x
eq1 = np.log10(Ni1)-np.log10(Ni0)-np.log10(NH3)-K1
eq2 = np.log10(Ni2)-np.log10(Ni1)-np.log10(NH3)-K2
eq3 = np.log10(Ni3)-np.log10(Ni2)-np.log10(NH3)-K3
eq4 = np.log10(Ni4)-np.log10(Ni3)-np.log10(NH3)-K4
eq5 = np.log10(Ni5)-np.log10(Ni4)-np.log10(NH3)-K5
eq6 = np.log10(Ni6)-np.log10(Ni5)-np.log10(NH3)-K6
eq7 = Ni0 + Ni1 + Ni2 + Ni3 + Ni4 + Ni5 + Ni6 - Ni_tot
eq8 = Ni1 + 2*Ni2 + 3*Ni3 + 4*Ni4 + 5*Ni5 + 6*Ni6 - NH3_tot[i]

return np.array([eq1, eq2, eq3, eq4, eq5, eq6, eq7, eq8])

Я использую цикл для расчета значений параметров при различных NH3_tot, что позволяет наблюдать изменение концентраций различных видов. Решение осуществляется с использованием следующего:

Код: Выделить всё

# Arrays for data storage
NH3_tot = np.arange(0, 6, 0.1)
Ni0_values = []
Ni1_values = []
Ni2_values = []
Ni3_values = []
Ni4_values = []
Ni5_values = []
Ni6_values = []
Solutions = []

# Speciation determination with least_square at various NH3_tot
for i in range(len(NH3_tot)):
x0 = [1,0.1,0.1,0.1,0.1,0.1,0.1,0.1]
result = least_squares(equations, x0, bounds=([0,0,0,0,0,0,0,0],[Ni_tot,Ni_tot,Ni_tot,Ni_tot,Ni_tot,Ni_tot,Ni_tot,Ni_tot]))

Solution = result.x
Solutions.append(Solution)
Ni0_value,Ni1_value,Ni2_value,Ni3_value,Ni4_value,Ni5_value,Ni6_value, NH3 = Solution

# Extract concentrations of species
Ni0_values.append(Ni0_value)
Ni1_values.append(Ni1_value)
Ni2_values.append(Ni2_value)
Ni3_values.append(Ni3_value)
Ni4_values.append(Ni4_value)
Ni5_values.append(Ni5_value)
Ni6_values.append(Ni6_value)

Проблема, с которой я столкнулся, заключается в том, что некоторые значения (обычно для низких и высоких значений NH3_tot) совершенно неверны, особенно не соблюдая уравнение. 7, так как сумма видов превышает Ni_tot.
Я пробовал перемещать границы, но, похоже, это не помогает, а крайние значения в целом неверны.
При печати результата всегда отображается:
сообщение: условие завершения xtol выполнено.
успех: True
Указывает на эффективность функции работает правильно.
Изменение химической системы приводит к аналогичным результатам.
Изменение записи уравнений сильно влияет на результат (если знаете, почему, пожалуйста, помогите !!), а использование np.log10 дает наиболее точный, но все равно неудовлетворительный результат.
Изменение метода разрешения функции наименьших_квадратов не улучшает результат.
Вот полный скрипт со всеми переменными:

Код: Выделить всё

import numpy as np
import matplotlib.pyplot as plt
from scipy.optimize import least_squares
import pandas as pd

Ni_tot = 1
K1 = 2.81
K2 = 2.27
K3 = 1.77
K4 = 1.27
K5 = 0.81
K6 = 0.15

# Define the system of nonlinear equations
def equations(x):
Ni0, Ni1, Ni2, Ni3, Ni4, Ni5, Ni6, NH3 = x
eq1 = np.log10(Ni1)-np.log10(Ni0)-np.log10(NH3)-K1
eq2 = np.log10(Ni2)-np.log10(Ni1)-np.log10(NH3)-K2
eq3 = np.log10(Ni3)-np.log10(Ni2)-np.log10(NH3)-K3
eq4 = np.log10(Ni4)-np.log10(Ni3)-np.log10(NH3)-K4
eq5 = np.log10(Ni5)-np.log10(Ni4)-np.log10(NH3)-K5
eq6 = np.log10(Ni6)-np.log10(Ni5)-np.log10(NH3)-K6
eq7 = Ni0 + Ni1 + Ni2 + Ni3 + Ni4 + Ni5 + Ni6 - Ni_tot
eq8 = Ni1 + 2*Ni2 + 3*Ni3 + 4*Ni4 + 5*Ni5 + 6*Ni6 - NH3_tot[i]

return np.array([eq1, eq2, eq3, eq4, eq5, eq6, eq7, eq8])

# Arrays for data storage
NH3_tot = np.arange(0, 14, 0.1)
Ni0_values = []
Ni1_values = []
Ni2_values = []
Ni3_values = []
Ni4_values = []
Ni5_values = []
Ni6_values = []
Solutions = []

# Speciation determination with least_square at various NH3_tot
for i in range(len(NH3_tot)):
x0 = [0.1,0.1,0.1,0.1,0.1,0.1,0.1,0.1]
result = least_squares(equations, x0, bounds=([0,0,0,0,0,0,0,0],[Ni_tot,Ni_tot,Ni_tot,Ni_tot,Ni_tot,Ni_tot,Ni_tot,Ni_tot]))
print(result)

Solution = result.x
Solutions.append(Solution)
Ni0_value,Ni1_value,Ni2_value,Ni3_value,Ni4_value,Ni5_value,Ni6_value, NH3 = Solution

# Extract concentrations of species
Ni0_values.append(Ni0_value)
Ni1_values.append(Ni1_value)
Ni2_values.append(Ni2_value)
Ni3_values.append(Ni3_value)
Ni4_values.append(Ni4_value)
Ni5_values.append(Ni5_value)
Ni6_values.append(Ni6_value)

# Conversion to DataFrame
results = pd.DataFrame(Solutions, columns=['Ni', 'Ni(NH3)', 'Ni(NH3)2', 'Ni(NH3)3', 'Ni(NH3)4', 'Ni(NH3)5', 'Ni(NH3)6', 'NH3'])

# Plot the concentrations
plt.figure(figsize=(10, 6))
plt.plot(NH3_tot, Ni0_values, label='Ni')
plt.plot(NH3_tot, Ni1_values, label='Ni1')
plt.plot(NH3_tot, Ni2_values, label='Ni2')
plt.plot(NH3_tot, Ni3_values, label='Ni3')
plt.plot(NH3_tot, Ni4_values, label='Ni4')
plt.plot(NH3_tot, Ni5_values, label='Ni5')
plt.plot(NH3_tot, Ni6_values, label='Ni6')

plt.xlabel('eq NH3')
plt.ylabel('Concentration')
plt.title('Concentration des espèces en fonction de NH3')
plt.legend()
plt.grid(True)
plt.show()

Цель — изучить вариацию распределения видов для различных лигандов, т.е. разных значений K. Но если видообразование сгенерировано неправильно, я не могу доверять своим данным.
Изменить после того, как помог jlandercy:
Используя исправление @jlandercy, нам удалось снизить неточность форму кривой с использованием математического выражения доли каждого комплекса (см. Alpha_i).
Однако, похоже, существует проблема, связанная с «концентрацией» лиганда. Действительно, физически не имеет смысла иметь зависимость образования комплекса (Ni-NH3) только от концентрации лиганда, а соотношение NH3/Ni_tot должно быть переменной величиной. Ожидаемый результат будет таким, как показано ниже:

С C_Ni_tot = 1 моль/л и C_NH3_tot = 0 -> 12 моль/л.
Используя сценарий @jlandercy со следующей модификацией:

Код: Выделить всё

L = np.linspace(0, 0.2, 100)
As = alphas(pK, L)

Мне удалось получить что-то того же типа, что показано:


Однако проблема «лигандного масштаба» все еще существует, что странно, и я не вижу этой проблемы.

Подробнее здесь: https://stackoverflow.com/questions/785 ... speciation