scipy 如何在python中从滴定点创建函数？[duplicate]

2admgd59  于 2022-11-10  发布在  Python
关注(0)|答案(1)|浏览(140)

此问题在此处已有答案：

Unable to fit curves to data points using curve_fit() from scipy because of "Optimal parameters not found" error [python]（1个答案）
三个月前关门了。
我希望在运行以下脚本时，通过在图形上指定点来自动计算函数。
滴定曲线与逻辑回归非常相似，但它们在转折点的斜率不同，因此不是100%相同的曲线。最好是有一个输出，它给出的函数与下面的函数类似：y = (A) / (1 + B^{x-C}) + D
顺便说一句，对不起，我的英语，这不是我的母语。请提出问题，如果它不清楚。


# from https://www.geeksforgeeks.org/graph-plotting-in-python-set-1/ inspired

# YOU HAVE TO HAVE INSTALLED THE FOLLOWING LIBRARIES : matplotlib.pyplot, colorama

# [pip install matplotlib.pyplot]

# [pip install colorama]

# normally alrady installed: math, time

import matplotlib.pyplot as plt
import math
from colorama import Fore
from colorama import Style
import time

### ANALYTE in ml ###

a = input(" Enter the constant volume of the analyte in ml: ")
print(f"{Fore.GREEN}The volume of the analyte is{Style.RESET_ALL}", a, f"{Fore.GREEN}ml{Style.RESET_ALL}")              
mol_a = input("Enter the molarity of the analyte: ")
print(f"{Fore.GREEN}The molarity of the analyte is{Style.RESET_ALL}", mol_a)

### TITRANT in ml ###

mol_t = input("Enter the molarity of the titrant: ")
t  = float(a)*float(mol_a) / float(mol_t)-0.00001
print(f"{Fore.GREEN}The molarity of the titrant is{Style.RESET_ALL}", mol_t)

# wait 1/2 a sec

print("calculating...")
time.sleep(.5)
print(f"{Fore.GREEN}The volume of the titrant under the equivalence point is{Style.RESET_ALL}",t, f"{Fore.GREEN}ml{Style.RESET_ALL}")

# Volumeof titrant in ml:

t0= 0.00001
t1= math.trunc(float(t)*10000)/10000
t2= float(t)*0.9999
t3= float(t)*0.999
t4= float(t)*0.99
t5= float(t)*0.97
t6= float(t)*0.9
t7= float(t)*0.7
t8= float(t)*0.5
t9= float(t)*0.3
t10= float(t)*0.1
t11= (float(t1))*1.0001
t12= (float(t1))*1.001
t13= (float(t1))*1.01
t14= (float(t1))*1.03
t15= (float(t1))*1.1
t16= (float(t1))*1.3
t17= (float(t1))*1.5
t18= (float(t1))*1.7
t19= (float(t1))*1.9
t20= (float(t1))*2

# mol of Analyte and Titrant:

anzahl_mol_von_starker_analyt = float(a)* 10**(-3)*float(mol_a)

# M                            =             z.B. 0.005      mol

anzahl_mol_von_starker_titrant = float(t)*10**(-3)*float(mol_t)

# pH at this ratio:

berechne_pH_t0       =     (-1)*math.log10(        (float(anzahl_mol_von_starker_analyt) -    float(float(t0)*10**(-3)*float(mol_t)))       /       (float(t0)*10**(-3)  + float(a)*10**(-3)   )       )
berechne_pH_t1       =     ( 7)
berechne_pH_unter_t2 =     (-1)*math.log10(        (float(anzahl_mol_von_starker_analyt) -    float(float(t2)*10**(-3)*float(mol_t)))       /       (float(t2)*10**(-3)  + float(a)*10**(-3)   )       )
berechne_pH_unter_t3 =     (-1)*math.log10(        (float(anzahl_mol_von_starker_analyt) -    float(float(t3)*10**(-3)*float(mol_t)))       /       (float(t3)*10**(-3)  + float(a)*10**(-3)   )       )
berechne_pH_unter_t4 =     (-1)*math.log10(        (float(anzahl_mol_von_starker_analyt) -    float(float(t4)*10**(-3)*float(mol_t)))       /       (float(t4)*10**(-3)  + float(a)*10**(-3)   )       )
berechne_pH_unter_t5 =     (-1)*math.log10(        (float(anzahl_mol_von_starker_analyt) -    float(float(t5)*10**(-3)*float(mol_t)))       /       (float(t5)*10**(-3)  + float(a)*10**(-3)   )       )
berechne_pH_unter_t6 =     (-1)*math.log10(        (float(anzahl_mol_von_starker_analyt) -    float(float(t6)*10**(-3)*float(mol_t)))       /       (float(t6)*10**(-3)  + float(a)*10**(-3)   )       )
berechne_pH_unter_t7 =     (-1)*math.log10(        (float(anzahl_mol_von_starker_analyt) -    float(float(t7)*10**(-3)*float(mol_t)))       /       (float(t7)*10**(-3)  + float(a)*10**(-3)   )       )
berechne_pH_unter_t8 =     (-1)*math.log10(        (float(anzahl_mol_von_starker_analyt) -    float(float(t8)*10**(-3)*float(mol_t)))       /       (float(t8)*10**(-3)  + float(a)*10**(-3)   )       )
berechne_pH_unter_t9 =     (-1)*math.log10(        (float(anzahl_mol_von_starker_analyt) -    float(float(t9)*10**(-3)*float(mol_t)))       /       (float(t9)*10**(-3)  + float(a)*10**(-3)   )       )
berechne_pH_unter_t10 =    (-1)*math.log10(        (float(anzahl_mol_von_starker_analyt) -    float(float(t10)*10**(-3)*float(mol_t)))      /       (float(t10)*10**(-3) + float(a)*10**(-3)   )       )
berechne_pH_ueber_t11 = 14-(-1)*math.log10(        (float(float(t11)*10**(-3)*float(mol_t) - (float(anzahl_mol_von_starker_analyt))))       /       (float(t11)*10**(-3) + float(a)*10**(-3)   )       )
berechne_pH_ueber_t12 = 14-(-1)*math.log10(        (float(float(t12)*10**(-3)*float(mol_t) - (float(anzahl_mol_von_starker_analyt))))       /       (float(t12)*10**(-3) + float(a)*10**(-3)   )       )
berechne_pH_ueber_t13 = 14-(-1)*math.log10(        (float(float(t13)*10**(-3)*float(mol_t) - (float(anzahl_mol_von_starker_analyt))))       /       (float(t13)*10**(-3) + float(a)*10**(-3)   )       )
berechne_pH_ueber_t14 = 14-(-1)*math.log10(        (float(float(t14)*10**(-3)*float(mol_t) - (float(anzahl_mol_von_starker_analyt))))       /       (float(t14)*10**(-3) + float(a)*10**(-3)   )       )
berechne_pH_ueber_t15 = 14-(-1)*math.log10(        (float(float(t15)*10**(-3)*float(mol_t) - (float(anzahl_mol_von_starker_analyt))))       /       (float(t15)*10**(-3) + float(a)*10**(-3)   )       )
berechne_pH_ueber_t16 = 14-(-1)*math.log10(        (float(float(t16)*10**(-3)*float(mol_t) - (float(anzahl_mol_von_starker_analyt))))       /       (float(t16)*10**(-3) + float(a)*10**(-3)   )       )
berechne_pH_ueber_t17 = 14-(-1)*math.log10(        (float(float(t17)*10**(-3)*float(mol_t) - (float(anzahl_mol_von_starker_analyt))))       /       (float(t17)*10**(-3) + float(a)*10**(-3)   )       )
berechne_pH_ueber_t18 = 14-(-1)*math.log10(        (float(float(t18)*10**(-3)*float(mol_t) - (float(anzahl_mol_von_starker_analyt))))       /       (float(t18)*10**(-3) + float(a)*10**(-3)   )       )
berechne_pH_ueber_t19 = 14-(-1)*math.log10(        (float(float(t19)*10**(-3)*float(mol_t) - (float(anzahl_mol_von_starker_analyt))))       /       (float(t19)*10**(-3) + float(a)*10**(-3)   )       )
berechne_pH_ueber_t20 = 14-(-1)*math.log10(        (float(float(t20)*10**(-3)*float(mol_t) - (float(anzahl_mol_von_starker_analyt))))       /       (float(t20)*10**(-3) + float(a)*10**(-3)   )       )

# X-axis: Volume of Titrant in ml:

x = [float(t0),float(t1),float(t2),float(t3),float(t4),float(t5),float(t6),float(t7),float(t8),float(t9),float(t10),float(t11),float(t12),float(t13),float(t14),float(t15),float(t16),float(t17),float(t18),float(t19),float(t20)]

# Y-axis: pH-value at this ratio:

y = [berechne_pH_t0,berechne_pH_t1,berechne_pH_unter_t2,berechne_pH_unter_t3,berechne_pH_unter_t4,berechne_pH_unter_t5,berechne_pH_unter_t6,berechne_pH_unter_t7,berechne_pH_unter_t8,berechne_pH_unter_t9,berechne_pH_unter_t10,berechne_pH_ueber_t11,berechne_pH_ueber_t12,berechne_pH_ueber_t13,berechne_pH_ueber_t14,berechne_pH_ueber_t15,berechne_pH_ueber_t16,berechne_pH_ueber_t17,berechne_pH_ueber_t18,berechne_pH_ueber_t19,berechne_pH_ueber_t20 ]

# Hiermit werden die Punkte in das Koordinatensystem geplottet:

plt.scatter(x, y, label= "values", color= "green", 
            marker= "o", s=30) 

# X-axis; title:

plt.xlabel('Volume of titrant in ml')

# Y-axis; title:

plt.ylabel('pH-value at this ratio')

# Title of graph:

plt.title('Titration curve, strong acid, strong base')

# show legend

plt.legend() 

# show plot

plt.show()
scipy

1条答案

按热度按时间
bq9c1y66

bq9c1y661#

为了使参数函数的参数与数据拟合，可以使用scipy.optimize.curve_fit：

def my_function(x, A, B, C, D, E):
    return A/(1 + B**(x - C)) + D + E*x

parameters, covariance = curve_fit(f = my_function, xdata = x, ydata = y)

如有必要，您可以通过传递参数的一些界限来改进结果：

parameters, covariance = curve_fit(f = my_function, xdata = x, ydata = y, bounds = ([0, 0, 0.9*t, -10, -10], [14, 1, 1.1*t, 10, 10]))

完整代码


# IMPORT

import matplotlib.pyplot as plt
import numpy as np
from scipy.optimize import curve_fit

# INPUT VALUES

a = 100
mol_a = 10
mol_t = 5

# DATA GENERATION

N = 40
t  = a*mol_a/mol_t

x1 = np.linspace(0, t - 1e-6, N//2)
x2 = np.linspace(t + 1e-6, 2*t, N//2)
anzahl_mol_von_starker_analyt = a*10**(-3)*mol_a
anzahl_mol_von_starker_titrant = t*10**(-3)*mol_t

y1 = (-1)*np.log10((anzahl_mol_von_starker_analyt - x1*10**(-3)*mol_t)/(x1*10**(-3) + a*10**(-3)))
y2 = 14 - (-1)*np.log10((x2*10**(-3)*mol_t - anzahl_mol_von_starker_analyt)/(x2*10**(-3) + a*10**(-3)))

x = np.concatenate((x1, x2))
y = np.concatenate((y1, y2))

# PARAMETERS FITTING

def my_function(x, A, B, C, D, E):
    return A/(1 + B**(x - C)) + D + E*x

parameters, covariance = curve_fit(f = my_function, xdata = x, ydata = y, bounds = ([0, 0, 0.9*t, -10, -10], [14, 1, 1.1*t, 10, 10]))

for parameter, name in zip(parameters, ['A', 'B', 'C', 'D', 'E']):
    print(f'{name} = {parameter:14.10f}')

x_fitted = np.linspace(x[0], x[-1], 1000)
y_fitted = my_function(x_fitted, *parameters)

# PLOT

plt.style.use('seaborn-darkgrid')
fig, ax = plt.subplots()

ax.plot(x, y, label = 'data', marker = 'o', linestyle = '')
ax.plot(x_fitted, y_fitted, label = 'fitted curve')

ax.set_xlabel('Volume of titrant in ml')
ax.set_ylabel('pH-value at this ratio')
ax.set_title('Titration curve, strong acid, strong base')
ax.legend(frameon = True)

plt.show()

结果

A =  13.2787117583
B =   0.7311691018
C = 199.8189983959
D =  -0.8940392051
E =   0.0054250609

图

赞(0）分享  回复(0）举报  2022-11-10
我来回答

相关问题

    查看更多

    热门标签

    更多
    JavaquerypythonNode开发语言requestUtil数据库Table后端算法LoggerMessageElementParser

    最新问答

    更多
    • 蜀ICP备13028337号-1 大数据知识库 https://www.saoniuhuo.com © All rights reserved
    • 本站内容来源互联网,如果侵犯您的权益请联系我们删除, 联系方式：448109455@qq.com