Python计算旋转圆盘电极LSV得到反应速率常数

会有学化学的看到吗？

只到拟合得到斜率和截距，后续计算很简单，暂时没写（懒得写）

书接上回：判断完LSV数据是否可用后，如果数据没问题就到计算速率常数这一步

主要需要获取的参数就是电流的-1次方和转速的-0.5次方线性拟合后得到的截距，只需要把截距存储起来，然后对截距取log，对过电位进行线性拟合即可。

Step1：输入转速，进行数学计算，存储输出

# 数据库导入
import pandas as pd
import numpy as np
import statsmodels.api as sm
import matplotlib.pyplot as plt
import math

# PART0:数据列表
potential_list = []#电压列表
rotating_speed_rpm=[]#转速列表rpm
rotating_speed_rad=[]#转速列表rad
rotating_speed_rad_half=[]#转速的-0.5次方列表
pl=potential_list
interception=[]#截距列表
interception_log=[]#截距列表
overpotentials=[]#过电位列表
R2list=[]#线性拟合数据方差

# PART1：输入相关参数
n= int(input("转速个数："))
m= int(input("取点数:"))
overpotential = input("过电位间隔(V):")# 输入过电位间隔
op= float(overpotential)

# PART2：输入转速rpm，输出为rad/s
r=1#循环
while r<=n:
    r=r+1
    rpm=int(input("输入转速rpm："))
    rad = rpm*0.10472
    rad_half=rad**-0.5
    rotating_speed_rpm.append(rpm)
    rotating_speed_rad.append(rad)
    rotating_speed_rad_half.append(rad_half)
print("转速",rotating_speed_rpm,"rpm")
print("转速",rotating_speed_rad,"rad/s")

Step2：数据进行查找计算存储输出，这里套用了一个循环语句，如果是对同一个电对，不同测试次数，对表格1进行完计算后，在对表格2进行计算，避免重复输入转速，如果两次数据转速不同则需要终止

p=int(input("开始输入1，结束输入其他数字:"))
while p==1:
    df = pd.read_excel(input("输入表格名称"))
    b = float(input("输入电压(V)：")) + op  # 输入平衡电压
    #PART3：循环检索表格中最接近的电压，并且输出所在行的数据，并对其进行数据处理
    a=1#循环计次用
    while a <= m:
    # 要搜索的数值

        search_column = 'Voltage'#确定检索的列
        search_value = b#检索的数字
    # 计算距离每个数值的差值，然后使用argmin找到最接近的索引
        closest_index = np.argmin(np.abs(df[search_column] - search_value))
    # 获取最接近的数据
        closest_value = df.loc[closest_index, search_column]#最接近的数据
        df7 = df[df['Voltage']==closest_value]  # 最接近数据
        print(df7)
        df8=1/df7#对电流取倒数
        arr = df8.iloc[:,1:n+1]#除去第一个数据
        arr_list1=arr.values.tolist()
        print("电流的倒数：",arr_list1[0])

        ops=-op*a*1000#过电位
        overpotentials.append(ops)

        

 

Step3：对数据进行线性拟合，输出相关数据

        y=arr_list1[0]
        x=rotating_speed_rad_half
        x_n = sm.add_constant(x)
        model = sm.OLS(y, x_n)
        coef = np.polyfit(x, y, 1)
        intercept = coef[1]
        intercept=float(intercept)
        interception.append(intercept)

        intercept_log = -math.log10(-intercept)
        interception_log.append(intercept_log)

        results=model.fit()
        R2=results.rsquared
        R2list.append(R2)
        # print('斜率', slope)
        print("截距", intercept)
        print("截距对数",intercept_log)
        print('R2: ', results.rsquared)

 Step4：输出得到的拟合结果，查看R方

        a = a+1#循环
        b = op+b#电压
        pl.append(closest_value)#记录电压



    print("所选电压:",pl)
    print("所有截距",interception)
    print("所有截距对数",interception_log)
    print("过电位",overpotentials,"mV")
    print("转速的-0.5次方和电流的-1次方线性拟合的方差",R2list)

 输出得到R方，可以看R方大小，R2太小则说明线性不好

 Step5：对log截距和过电位线性拟合

    y=overpotentials
    x=interception_log
    excel = []
    excel.append(x)
    excel.append(y)
    df = pd.DataFrame(excel)
    df2 = df.transpose()
    x_n = sm.add_constant(x)  # statsmodels进行回归时，一定要添加此常数项
    coefficent = np.polyfit(x, y, 1)
    y_fixed=np.poly1d(coefficent)
    model = sm.OLS(y, x_n)  # sm.OLS是回归分析模型
    results = model.fit()  # results是回归分析后的结果
    print(results.summary())
    print('Parameters: ', results.params)
    print('R2: ', results.rsquared)

 Step6:出图

    plt.rcParams['font.sans-serif'] = [u'simHei']  # 显示中文
    plt.rcParams['axes.unicode_minus'] = False  # 解决负号问题
    plt.rcParams['font.size'] = 20
    plt.rcParams['font.family'] = 'Arial'
    plt.figure(figsize=(8, 6))  # 单位为英寸

    plt.title(u"KL-Curves",fontsize=20)
    plt.xlabel(u"-logik",fontsize=20)
    plt.ylabel(u"overpotential(mV)",fontsize=20)



    plt.scatter(x, y, marker="o",color="b", s=50)
    plt.plot(x, y_fixed(x), linewidth=3, color="r")
    plt.show()

 斜率截距

 

Step7：对表格进行清空，便于下一次计算，对应Step2

    potential_list.clear()
    interception.clear()
    interception_log.clear()
    overpotentials.clear()
    arr_list1.clear()
    R2list = []  # 线性拟合数据方差

Step8：是否输出文件

    output = int(input("是否输出结果，是输入1，否0："))
    if output == 1:

        df2.to_excel(input("输出：输入文件名："), index=False)
    p=int(input("继续输入1，结束输入其他数字"))
    a=0