Python日记（2）——反应速率的计算

Python日记（2）——反应速率的计算

每天做一个Python小练习，顺便记录一些小技巧。
化工数值计算实验的实验二：1，3-丁二烯气相二聚反应的速率计算
一、目的和要求
1.掌握常用的插值方法
2.掌握数值微分的计算方法
二、问题提出
1.函数插值
函数是工程实践中重要的工具，它能够描述变量之间的关系。但是实际应用
当中，我们很难使用解析表达式描述两变量之间的关系。我们能够得到的，只是
从实验测量得到的离散的数据，即函数表。当我们需要使用的数据点不在函数表
中时，我们就需要进行插值。
2.数值微分
有时候，需要从实验得到的有限的数据中推断出变量之间的相对变化率。但
是，仅从有限个数据，不可能通过解析的方法精确求解出其导数，只能利用数值
方法计算。数值微分与插值和拟合是密不可分的。在进行数值计算时，往往需要
插值或拟合来增加数据点的密度，提高精度。
3.二聚反应速率的计算

t/min	p/mmHg
0	632.0
5	611.0
10	592.0
15	573.5
20	558.5
25	545.0
30	533.0
35	523.0
40	514.0
45	505.0
50	497.0
55	490.0
60	484.0
65	478.5
70	473.0
75	468.0
80	463.0
90	458.0
95	453.0

丁二烯可发生二聚反应：
2C4H6 → (C4H6)2
若初始压力为p0，在t时刻丁二烯的分压为p1，二聚体的分压为p2，则它们存
在关系p1+ 2p2= p0。该反应的速率可表示为
-1/2dp1/dt = dp2/dt
也就是“虽计+算反应速率，需要知道P或p2与时间的函数关系，但是实验只能测量总压，而且实验测量只能得到有限的离散的数据，因此需要从有限的数据中计算分压对时间的导数，从而得到反应速率。
三、实验内容
1.换算p的单位
根据 1kPa = 7.5mmHg

def p_convert(p):
    for i in p:
        index_num = p.index(i)
        i /= 7.5
        p[index_num] = round(i, 4)
    return p

2.因为二聚体的分压等于Δp，所以利用循环求Δp并返回便是p1的分压

def p_delta(p):
    p1 = list()
    for i in range(0, len(p)):
        if i+1 == len(p):
            break
        else:
            pl = round(p[0] - p[i])
            p1.append(pl)
    return p1

3.曲线插值
（1）设置x的点集
linspace(begin, end, num)
begin为左端点，end为右端点，num为元素个数

def generate_data(begin, end, num):
    """ 产生x点集 """
    x = np.linspace(begin, end, num)
    return x

（1）线性插值
interp1d(x, y, kind)
x, y为已知数据点，kind为可选类型，可以是字符串或整数

def interpolate_linear(x, y):
    """ 线性插值 """
    f_linear = interpolate.interp1d(x, y)
    return f_linear

(3)求得样条曲线或导数，默认为0阶导数。
scipy.interpolate.splrep(x, y, w=None, xb=None, xe=None, k=3, task=0, s=None, t=None, full_output=0, per=0, quiet=1)
查找一维曲线的B样条曲线表示。给定数据点集，确定区间上度k的平滑样条近似。
（各个参数的表达含义我还没学会）
x, y为已知数据点，x_new就是新的细分后的x点集，元素个数要与差值后的个数保持一直。

def interpolate_b