CN110514608B - 一种基于光谱的反应动力学速率常数的无偏估计方法 - Google Patents

一种基于光谱的反应动力学速率常数的无偏估计方法 Download PDF

Info

Publication number
CN110514608B
CN110514608B CN201910802085.8A CN201910802085A CN110514608B CN 110514608 B CN110514608 B CN 110514608B CN 201910802085 A CN201910802085 A CN 201910802085A CN 110514608 B CN110514608 B CN 110514608B
Authority
CN
China
Prior art keywords
matrix
rate constant
reaction
spectrum
reaction kinetic
Prior art date
Legal status (The legal status is an assumption and is not a legal conclusion. Google has not performed a legal analysis and makes no representation as to the accuracy of the status listed.)
Active
Application number
CN201910802085.8A
Other languages
English (en)
Other versions
CN110514608A (zh
Inventor
陈伟锋
何腾
张贵军
胡俊
Current Assignee (The listed assignees may be inaccurate. Google has not performed a legal analysis and makes no representation or warranty as to the accuracy of the list.)
Zhejiang University of Technology ZJUT
Original Assignee
Zhejiang University of Technology ZJUT
Priority date (The priority date is an assumption and is not a legal conclusion. Google has not performed a legal analysis and makes no representation as to the accuracy of the date listed.)
Filing date
Publication date
Application filed by Zhejiang University of Technology ZJUT filed Critical Zhejiang University of Technology ZJUT
Priority to CN201910802085.8A priority Critical patent/CN110514608B/zh
Publication of CN110514608A publication Critical patent/CN110514608A/zh
Application granted granted Critical
Publication of CN110514608B publication Critical patent/CN110514608B/zh
Active legal-status Critical Current
Anticipated expiration legal-status Critical

Links

Images

Classifications

    • G—PHYSICS
    • G01—MEASURING; TESTING
    • G01N—INVESTIGATING OR ANALYSING MATERIALS BY DETERMINING THEIR CHEMICAL OR PHYSICAL PROPERTIES
    • G01N21/00—Investigating or analysing materials by the use of optical means, i.e. using sub-millimetre waves, infrared, visible or ultraviolet light
    • G01N21/17—Systems in which incident light is modified in accordance with the properties of the material investigated
    • G01N21/25—Colour; Spectral properties, i.e. comparison of effect of material on the light at two or more different wavelengths or wavelength bands
    • G01N21/27—Colour; Spectral properties, i.e. comparison of effect of material on the light at two or more different wavelengths or wavelength bands using photo-electric detection ; circuits for computing concentration
    • G01N21/274—Calibration, base line adjustment, drift correction
    • G—PHYSICS
    • G01—MEASURING; TESTING
    • G01N—INVESTIGATING OR ANALYSING MATERIALS BY DETERMINING THEIR CHEMICAL OR PHYSICAL PROPERTIES
    • G01N21/00—Investigating or analysing materials by the use of optical means, i.e. using sub-millimetre waves, infrared, visible or ultraviolet light
    • G01N21/17—Systems in which incident light is modified in accordance with the properties of the material investigated
    • G01N21/25—Colour; Spectral properties, i.e. comparison of effect of material on the light at two or more different wavelengths or wavelength bands
    • G01N21/31—Investigating relative effect of material at wavelengths characteristic of specific elements or molecules, e.g. atomic absorption spectrometry
    • G—PHYSICS
    • G01—MEASURING; TESTING
    • G01N—INVESTIGATING OR ANALYSING MATERIALS BY DETERMINING THEIR CHEMICAL OR PHYSICAL PROPERTIES
    • G01N21/00—Investigating or analysing materials by the use of optical means, i.e. using sub-millimetre waves, infrared, visible or ultraviolet light
    • G01N21/84—Systems specially adapted for particular applications
    • G01N2021/8411—Application to online plant, process monitoring

Landscapes

  • Physics & Mathematics (AREA)
  • Spectroscopy & Molecular Physics (AREA)
  • General Physics & Mathematics (AREA)
  • Immunology (AREA)
  • Chemical & Material Sciences (AREA)
  • Analytical Chemistry (AREA)
  • Biochemistry (AREA)
  • General Health & Medical Sciences (AREA)
  • Health & Medical Sciences (AREA)
  • Life Sciences & Earth Sciences (AREA)
  • Pathology (AREA)
  • Engineering & Computer Science (AREA)
  • Mathematical Physics (AREA)
  • Theoretical Computer Science (AREA)
  • Investigating Or Analysing Materials By Optical Means (AREA)
  • Investigating Or Analysing Materials By The Use Of Chemical Reactions (AREA)

Abstract

本发明提供了一种基于光谱的反应动力学速率常数的无偏估计方法，属于化学反应动力学领域。该方法首先构建光谱数据结构，根据反应物初始浓度和反应动力学速率常数的初始估计值，计算浓度分布曲线，获得将光谱矩阵，与测量得到的光谱矩阵比较，得到误差矩阵；然后根据误差矩阵调整反应动力学速率常数的估计值，重新计算误差矩阵，重复迭代，直到误差达到最小值；最后用最终误差矩阵估计噪声方差，在测量得到的光谱中增加噪声，重新计算反应动力学速率常数估计值，使用无偏估计方法，得到的结果更加接近真实值。由于使用无偏估计方法，不依赖反应物质的纯光谱等其他信息，有效的减少了噪声产生的影响，使得到的反应动力学速率常数更加精确。

Description

一种基于光谱的反应动力学速率常数的无偏估计方法
技术领域
本发明属于化学反应动力学领域，具体地涉及一种基于光谱的反应动力学速率常数的无偏估计方法。
背景技术
建立精确的化学反应过程的动力学模型，以便对化学反应过程进行合理的预测和有效的控制，既能保证整个生产过程的安全性，也能提高化工生产效率，降低成本，最大化经济效益，对于现代化学工业具有重大意义。精确的反应动力学模型包括了模型结构和模型参数。确立模型结构不仅需要对化学反应过程的反应机理十分了解，还需要根据反应器的类型，所研究问题的特点进行合理的假设。同样十分重要的是如何优化模型参数，使得模型模拟得出的数据与实际值偏差最小，是当前化学反应动力学研究的一个重点。
经过多年的发展，大量学者从优化设定及最优控制等方面对模型参数估计展开了研究。基于光谱拟合动力学模型参数的方法也不断发展，并针对影响结果精度的因素提出了相应的策略。经典曲线分析方法是一种非线性最小二乘拟合方法，该方法是一种梯度方法，从初始值出发依据梯度信息迭代到最优值。使用该方法可以较为快速地从光谱中估计出反应速率常数，但是没有对实验误差和测量噪声进行进一步的研究。为了补充该缺陷，利用光谱估算反应速率常数的约束最小二乘法也被提出，但是该方法在参与反应物质的纯光谱难以获得的时候，效果大打折扣。因此现有的反应动力学速率常数估计方法仍存在缺陷，需要改进。
发明内容
为了解决现有的反应动力学参数估计方法存在的不足，本发明提出一种基于光谱的反应动力学速率参数的无偏估计方法。
本发明解决其技术问题所采用的技术方案是：一种基于光谱的反应动力学速率常数的无偏估计方法，具体包括以下步骤：
步骤一：基于Beer-Lambert定律构建光谱数据的结构，使用常微分方程组建立反应动力学模型，根据反应物的初始浓度和反应动力学速率常数初始估计值，计算化学反应中涉及的所有物质随时间变化的浓度分布曲线，利用最小二乘计算各物质的吸光系数，计算光谱数据矩阵，将计算得到的光谱数据矩阵与实验测量的光谱比较，得到误差矩阵E(k0)。
步骤二：根据步骤一所述误差矩阵E(k0)调整反应动力学速率常数初始估计值，计算得到新的浓度分布曲线，利用最小二乘计算各物质的吸光系数，计算新的光谱数据矩阵，将计算得到的新光谱数据矩阵与实验测量的光谱比较，重新计算得到新的误差矩阵E(k1)。经上述方法重复迭代，获得最终误差矩阵E(kn-1)，直到最终误差矩阵E(kn-1)中每个元素的平方和与上一次误差矩阵E(kn-2)中每个元素的平方和之差小于10-6。n为迭代次数，n≥1。
步骤三：根据步骤二得到的最终误差矩阵E(kn-1)中每个元素的平方和估计噪声的方差，在所述测量得到的光谱中增加具有所述方差的噪声，获得具有叠加噪声的光谱数据，重新计算反应动力学速率常数估计值，使用无偏估计方法，完成对反应动力学速率常数的估计。
进一步地，步骤一具体为：
所述光谱数据矩阵D为：D＝CST+E(kn-1)，其中，C表示具有光吸收特性的物质随时间变化的浓度矩阵，S是每种物质在不同测量波长下的吸光系数矩阵，D表示反应溶液对于每种测量波长的吸光程度随时间的变化，E(kn-1)与D维度相同。
以时间为自变量的常微分方程组建立反应动力学模型，给定光谱数据矩阵D和所有反应物的初始浓度，根据反应动力学速率常数的初始估计值k0，计算出C、S、E(kn-1)，其中，C+为C的伪逆：
E(kn-1)＝D-CST
ST＝C+D
C+＝(CTC)-1CT
进一步地，步骤二包括以下子步骤：
(2.1)E(k0)中每个元素Ei,j的平方和ssq的求解过程为：
Figure BDA0002182608480000021
nt是采样点的数量，nλ为测量波长的数量。
(2.2)以ssq向其最小值移动的方式计算反应动力学速率常数k的偏移Δk，反应动力学速率常数初始值表示为k0，经过Δk＝(Δk1,Δk2,…,Δknk)的微小偏移后，用泰勒展开来估计E(k0)，其中E(k0)的泰勒展开式中含有偏导数矩阵。
(2.3)将E(k0)和所述偏导数矩阵矢量化，其中，矢量化后的偏导数矩阵称为雅可比矩阵J，通过下列公式求解Δk，e(k0)为E(k0)矢量化后的向量：
e(k0)＝-JΔk
J+＝(JTJ)-1JT
Δk＝-J+e(k0)
(2.4)将计算所得的Δk与k0相加，得到新的反应动力学速率常数估计值，重复步骤(2.1)-(2.3)，直到E(kn-1)的ssq和与E(kn-2)的ssq之差小于10-6，循环终止，获得此时的反应动力学速率常数k(D)。
进一步地，所述步骤三包括以下子步骤：
(3.1)计算步骤(2.4)中得到的速率常数k(D)所对应的E(kn-1)含有的所有元素的方差σ2，视为实验测量噪声的方差。
(3.2)在测量光谱数据矩阵D上加上均值为0、方差为σ2的高斯白噪声，记为Z，重新计算反应动力学速率常数的估计值k(Z)，得到无偏估计结果为k＝2k(D)-k(Z)。
(3.3)重复步骤(3.2)，将结果取平均，得到最终反应动力学速率常数无偏估计结果。
本发明的有益效果是：本发明从光谱出发，使用经典曲线分析方法实现对反应动力学速率常数的估计。由于实验中存在测量噪声，得到的结果含有一定的误差，因此使用无偏估计方法，不依赖反应物质的纯光谱等其他信息，有效的减少了噪声产生的影响，使得到的反应动力学速率常数更加精确。
附图说明
图1为本发明方法的流程示意图；
图2为步骤(2.2)中有限差分自适应步长选择方法的流程示意图。
具体实施方式
为了使本发明的目的和技术方案更加清楚明白，以下结合附图和实施例,对本发明进行进一步详细说明。应当理解，此处所描述的具体实施例仅仅用以解释本发明，并不用于限定本发明。
参照图1，一种基于光谱的反应动力学速率常数的无偏估计方法，包括以下步骤：
步骤一：基于Beer-Lambert定律构建光谱数据的结构，使用常微分方程组建立反应动力学模型，根据反应物的初始浓度和反应动力学速率常数初始估计值，计算化学反应中涉及的所有物质随时间变化的浓度分布曲线，计算化学反应中涉及的所有物质随时间变化的浓度分布曲线，利用最小二乘计算各物质的吸光系数，得到光谱数据矩阵，将计算得到的光谱矩阵与实验测量的光谱比较，得到误差矩阵E(k0)。
步骤二：根据步骤一所述误差矩阵调整反应动力学速率常数初始估计值，计算得到新的浓度分布曲线，利用最小二乘计算各物质的吸光系数，得到新光谱数据矩阵，将计算得到的新光谱数据矩阵与所述测量得到的光谱比较，重新计算得到新的误差矩阵E(k1)。经上述方法重复迭代，获得最终的误差矩阵E(kn-1)，直到最终的误差矩阵E(kn-1)中每个元素的平方和与上一次误差矩阵E(kn-2)中每个元素的平方和之差小于10-6。n为迭代次数，n≥1。
步骤三：根据步骤二得到的E(kn-1)中每个元素的平方和估计噪声的方差，在所述实验测量的光谱中增加具有所述方差的噪声，获得具有叠加噪声的光谱数据，重新计算反应动力学速率常数估计值，使用无偏估计方法，完成对反应动力学速率常数的估计。
其中，步骤一的方法具体为：根据Beer-Lambert定律，时间为t时，在特定的波长λ处，溶液的吸光度是来自所有吸收该波长光线的物质的吸光总和。每种物质对吸光度的贡献与浓度和光程长度(入射光束在溶液中行进的距离)成线性比例。如果路径长度和比例常数由单个常数si(λ)表示，则Beer-Lambert定律由以下公式表示：
d(t,λ)＝c1(t)s1(λ)+c2(t)s2(λ)+…+cnc(t)snc(λ)
nc为在波长λ处吸光物质的数量，ci(t)是时间为t时i物质的浓度，d(t,λ)表示时间为t时，溶液对波长为λ的光的吸收总和。跟踪化学反应的变化，获得的光谱即为反应时间的函数，测量结果为二维光谱数据矩阵：D＝CST+E(kn-1)，其中，C表示每种物质随时间变化的浓度，为nt×nc的矩阵；S是每种物质在不同测量波长下的吸光系数，为nλ×nc的矩阵；D为nt×nλ光谱数据矩阵，表示反应溶液对于每种测量波长的吸光程度随时间的变化；由于测量存在误差，引入E(kn-1)，与D维度相同。
本发明研究的化学反应，假定反应器内各种物料均匀混合，该假设下研究对象的模型用以时间为自变量的常微分方程组(ordinary differential equations，ODES)建立反应动力学模型，给定光谱数据矩阵D和所有参与反应物质的初始浓度，根据反应动力学速率常数的初始估计值k0，结合ODES，计算出C、E(kn-1)，其中，C+为C的伪逆：
E(kn-1)＝D-CST
ST＝C+D
C+＝(CTC)-1CT
步骤二还包括以下子步骤：
(2.1)使用经典曲线分析方法，计算E(k0)中每个元素Ei,j的平方和ssq(sum ofsquares，ssq)的求解过程为：
Figure BDA0002182608480000051
nt是采样点的数量，nλ为测量波长的数量。
(2.2)以ssq向其最小值移动的方式计算反应动力学速率常数k的偏移Δk，为此，需要强调的是E(k0)和随后的ssq仅是反应动力学速率常数k的函数，反应动力学速率常数初始值表示为k0，经过Δk＝(Δk1,Δk2,…,Δknk)的微小偏移后，用泰勒展开来估计E(k0)，其中E(k0)的泰勒展开式中含有偏导数矩阵，仅使用一阶导数，得到如下线性表达式：
Figure BDA0002182608480000052
虽然上述公式为近似表示方法，由于其为线性表达式，这使得它易于处理。目标是使E(k0+Δk)的值趋于0，因此用0代替E(k0+Δk)，上式重新排列结果如下：
Figure BDA0002182608480000053
误差矩阵E(k0)可以由上式表达，并且偏导数可以根据下式方法计算：
Figure BDA0002182608480000061
这里使用的是前向差分，本文根据Gill提出的算法进行选择Δki，算法流程参照图2，步骤详解如下：
Figure BDA0002182608480000062
Figure BDA0002182608480000063
Figure BDA0002182608480000064
Figure BDA0002182608480000065
Figure BDA0002182608480000066
Figure BDA0002182608480000067
εA＝|f(x)|εR
其中，h为有限差分自适应步长，
Figure BDA0002182608480000068
为函数f(x)的前向差分，
Figure BDA0002182608480000069
为函数f(x)的后向差分，
Figure BDA00021826084800000610
为函数f(x)的中心差分，Φ(f,h)为函数f(x)的二阶导数，εR为机器精度，εA为最小步长，
Figure BDA00021826084800000611
为第一判别变量，
Figure BDA00021826084800000612
为第二判别变量，
Figure BDA00021826084800000613
为第三判别变量。
步骤1：确认函数f(x)与自变量x，选择参数M，令m＝0，hs＝-1，
Figure BDA00021826084800000614
Figure BDA00021826084800000615
如下式：
Figure BDA00021826084800000616
计算式(9)中的各个参数；如果
Figure BDA00021826084800000617
则hs←h0。
步骤2：m＝m+1，hm＝10hm-1，重新计算式(9)中的各个参数。如果hs<0并且
Figure BDA00021826084800000618
则hs←hm。如果满足
Figure BDA00021826084800000619
则hΦ←hm，步骤2停止，进入下一步；否则重复步骤2直到m＝M，进入下一步；
步骤3：m＝m+1，hm＝hm-1/10，重新计算式(9)中的各个参数。如果
Figure BDA00021826084800000620
则hs←hm。如果
Figure BDA00021826084800000621
则hΦ←hm，步骤3停止，进入下一步；否则重复步骤3直到m＝M，步骤3停止，进入下一步；
步骤4：hF＝hΦ，输出hF，算法终止；
步骤5：如果hs<0，则取
Figure BDA0002182608480000071
如果
Figure BDA0002182608480000072
则取hF←hs；如果
Figure BDA0002182608480000073
则取hF←hM。输出hF，算法终止。
(2.3)将E(k0)和所述偏导数矩阵矢量化，即展开成长列向量，将所述长列向量折叠成矩阵与向量的乘积；矢量化后的偏导数矩阵称为雅可比矩阵J，通过下列公式求解Δk，e(k0)为E(k0)矢量化后的向量：
e(k0)＝-JΔk
J+＝(JTJ)-1JT
Δk＝-J+e(k0)
(2.4)将计算所得的Δk与k0相加，得到新的反应动力学速率常数估计值，重复步骤(2.1)-(2.3)，直到E(kn-1)的ssq和与E(kn-2)的ssq之差小于10-6，循环终止，获得此时的反应动力学速率常数k(D)。
步骤三还包括以下子步骤：
(3.1)计算步骤(2.4)中得到的速率常数k(D)所对应的E(kn-1)含有的所有元素的方差σ2，视为实验测量噪声的方差。
(3.2)在测量光谱数据矩阵D上加上均值为0、方差为σ2的高斯白噪声，记为Z，重新计算反应动力学速率常数的估计值k(Z)，得到无偏估计结果为k＝2k(D)-k(Z)。
(3.3)需要强调的是，无偏性的实际意义是指没有***性的偏差。无论测量几次同一化学反应过程，经过参数估计得到的反应速率常数相对于真实值而言，总有一些偏高，一些偏低。而无偏性表示，把这些误差平均起来，其值为零，即无偏估计量只存在随机误差。因此重复步骤(3.2)150次，将结果取平均，得到最终反应动力学速率常数无偏估计结果。
实施例
(1)根据Beer-Lambert定律，时间为t时，在特定的波长λ处，溶液的吸光度是来自所有吸收该波长光线的物质的吸光总和。每种物质对吸光度的贡献与浓度和光程长度(入射光束在溶液中行进的距离)成线性比例。跟踪化学反应的变化，获得的光谱即为反应时间的函数，测量结果是一个二维矩阵，即为光谱数据矩阵：D＝CST+E(kn-1)，其中，C表示每种物质随时间变化的浓度，S是每种物质在不同测量波长下的吸光系数，D为光谱数据矩阵，表示反应溶液对于每种测量波长的吸光程度随时间的变化，E(kn-1)表示误差矩阵，与D维度相同。
以时间为自变量的ODES建立反应动力学模型，这里选用MCR-ALS工具包中的演示算例，其ODES如下：
Figure BDA0002182608480000081
Figure BDA0002182608480000082
反应速率常数真实值kr＝[2,0.2]，物质初始浓度[A]0＝0.001mol·L-1，其余物质为0。采样时间t＝10s，t0＝0，采样点个数nt＝300。测量光谱矩阵D已知，反应动力学参数的初始估计值设为k0＝[1,0.5]，结合ODES可以计算出三种物质浓度随时间分布的矩阵C，误差矩阵E的计算方法如下，其中C+为C的伪逆：
E＝D-CST
ST＝C+D
C+＝(CTC)-1CT
(2.1)使用经典曲线分析方法，计算所述误差矩阵E(k0)中每个元素Ei,j的平方和ssq的求解过程为：
Figure BDA0002182608480000083
nt是采样点的数量，nλ为测量波长的数量。
(2.2)一旦确定了ssq，就以ssq向其最小值移动的方式计算反应速率常数的偏移。为此，需要强调的是E(k0)和随后的ssq仅是速率常数k的函数，如果速率常数初始值表示为k0，则参数经过Δk＝(Δk1,Δk2,…Δknk)的微小偏移之后，可以用泰勒展开来估计E(k0)，如果仅使用一阶导数，得到如下线性表达式：
Figure BDA0002182608480000084
其中
Figure BDA0002182608480000085
为E(k0)对于knk的偏导数矩阵。目标是使E(k0+Δk)的值趋于0，因此用0代替E(k0+Δk)，上式重新排列结果如下：
Figure BDA0002182608480000086
误差矩阵E(k0)可以由上式表达，并且偏导数可以根据下式方法计算：
Figure BDA0002182608480000087
(2.3)将E(k0)和所述偏导数矩阵矢量化，即展开成长列向量，将所述长列向量折叠成矩阵与向量的乘积；矢量化后的偏导数矩阵称为雅可比矩阵J，通过下列公式求解Δk，e(k0)为E(k0)矢量化后的向量：
e(k0)＝-JΔk
J+＝(JTJ)-1JT
Δk＝-J+e(k0)
(2.4)将计算所得的Δk与k0相加，得到新的反应动力学速率常数估计值，重复步骤(2.1)-(2.3)，直到E(kn-1)的ssq和与E(kn-2)的ssq之差小于10-6，循环终止，获得此时的反应动力学速率常数。
通常，如上所述方法在最小值附近快速收敛。然而，如果初始估计较差，则泰勒级数展开的函数逼近和问题的线性化可能导致ssq的发散和算法的失败。Marquardt基于Levenberg的思想提出的修改是在参数计算过程中将一定数量的Marquardt参数mp加到Hessian矩阵的对角元素上。Hessian矩阵及Δk新的表达式如下：
H＝JTJ
Δk＝-(H+mp*I)-1JTe(k0)
当mp的值显著大于Hessian矩阵的对角元素时，Δk表达式如下，具有线性下降的性质：
Figure BDA0002182608480000091
当mp很小时，Δk表达式不变。管理参数mp策略如下：最初设置为1，如果参数ssq开始发散，则mp不断增加(每次迭代乘以10)，ssq收敛则不断减小(每次迭代除以√10)，在迭代结束后重置为1。
(3.1)在步骤(2.4)中得到的速率常数记为k(D)＝[1.9766,0.19476]，计算它所对应的E(kn-1)含有的所有元素的方差σ2＝4.3，可以视为实验测量噪声的方差；
(3.2)在测量光谱矩阵D上加上均值为0，方差为σ2的高斯白噪声，记为Z，重新计算反应动力学速率常数的估计值k(Z)，得到无偏估计结果为k＝2k(D)-k(Z)。
(3.3)重复步骤(3.2)150次，结果取平均，算作一次无偏估计计算。以真实值kr＝[2,0.2]为标准，平均误差率可以作为评价一个结果是否优秀的指标。
未考虑测量噪声所得的结果为k(D)＝[1.9766,0.19476]，平均误差率为1.90％。使用了无偏估计方法之后，在100次的重复实验中，有69次k(Z)比k(D)更接近真实值，部分结果如表1。说明在此算例中，无偏估计能够较为有效地减少噪声带来的误差。
表1 MCR-ALS算例仿真无偏估计部分结果
Figure BDA0002182608480000101
以上显示和描述了本发明的基本原理和主要特征和本发明的优点，但本发明并不仅仅限于说明书和实施方式所列运用，它完全可以被适用于各种适合本发明的领域，对于熟悉本领域的人员而言，可容易地实现另外的修改，因此在不背离权利要求及等同范围所限制的一般概念下，本发明并不限于特定的细节和这里示出的实施例。

Claims (2)

1.一种基于光谱的反应动力学速率常数的无偏估计方法，其特征在于：具体包括以下步骤：
步骤一：基于Beer-Lambert定律构建光谱数据的结构，使用常微分方程组建立反应动力学模型，根据反应物的初始浓度和反应动力学速率常数初始估计值，计算化学反应中涉及的所有物质随时间变化的浓度分布曲线，利用最小二乘计算各物质的吸光系数，计算光谱数据矩阵，将计算得到的光谱数据矩阵与实验测量的光谱比较，得到误差矩阵
Figure DEST_PATH_IMAGE001
；
步骤二：根据步骤一所述误差矩阵
Figure 760067DEST_PATH_IMAGE001
调整反应动力学速率常数初始估计值，计算得到新的浓度分布曲线，利用最小二乘计算各物质的吸光系数，计算新的光谱数据矩阵，将计算得到的新光谱数据矩阵与实验测量的光谱比较，重新计算得到新的误差矩阵
Figure DEST_PATH_IMAGE002
；经上述方法重复迭代，获得最终误差矩阵
Figure DEST_PATH_IMAGE003
，直到最终误差矩阵
Figure 631202DEST_PATH_IMAGE003
中每个元素的平方和与上一次误差矩阵
Figure DEST_PATH_IMAGE004
中每个元素的平方和之差小于10-6；n为迭代次数，n≥1；
步骤二包括以下子步骤：
（2.1）
Figure 269994DEST_PATH_IMAGE001
中每个元素
Figure DEST_PATH_IMAGE005
的平方和ssq的求解过程为：
Figure DEST_PATH_IMAGE006
，
Figure DEST_PATH_IMAGE007
是采样点的数量，
Figure DEST_PATH_IMAGE008
为测量波长的数量；
（2.2）以ssq向其最小值移动的方式计算反应动力学速率常数
Figure DEST_PATH_IMAGE009
的偏移
Figure DEST_PATH_IMAGE010
，反应动力学速率常数初始值表示为
Figure DEST_PATH_IMAGE011
，经过
Figure DEST_PATH_IMAGE012
的微小偏移后，用泰勒展开来估计
Figure 163738DEST_PATH_IMAGE001
，其中
Figure 509269DEST_PATH_IMAGE001
的泰勒展开式中含有偏导数矩阵；
（2.3）将所述
Figure 66152DEST_PATH_IMAGE001
和所述偏导数矩阵矢量化，其中，矢量化后的偏导数矩阵称为雅可比矩阵
Figure DEST_PATH_IMAGE013
，通过下列公式求解
Figure 696635DEST_PATH_IMAGE010
，
Figure DEST_PATH_IMAGE014
为
Figure 364377DEST_PATH_IMAGE001
矢量化后的向量：
Figure DEST_PATH_IMAGE015
（2.4）将计算所得的
Figure 220206DEST_PATH_IMAGE010
与
Figure 885674DEST_PATH_IMAGE011
相加，得到新的反应动力学速率常数估计值，重复步骤（2.1）-（2.3），直到
Figure 797260DEST_PATH_IMAGE003
的ssp和与
Figure 65430DEST_PATH_IMAGE004
的ssp之差小于10-6，循环终止，获得此时的反应动力学速率常数
Figure DEST_PATH_IMAGE016
；
步骤三：根据步骤二得到的最终误差矩阵
Figure 979029DEST_PATH_IMAGE003
中每个元素的平方和估计噪声的方差，在所述测量得到的光谱中增加具有所述方差的噪声，获得具有叠加噪声的光谱数据，重新计算反应动力学速率常数估计值，使用无偏估计方法，完成对反应动力学速率常数的估计；
步骤三包括以下子步骤：
（3.1）计算步骤（2.4）中得到的速率常数
Figure 815398DEST_PATH_IMAGE016
所对应的
Figure 588182DEST_PATH_IMAGE003
含有的所有元素的方差
Figure DEST_PATH_IMAGE017
，视为实验测量噪声的方差；
（3.2）在测量光谱数据矩阵上加上均值为0、方差为
Figure 548791DEST_PATH_IMAGE017
的高斯白噪声，记为Z，重新计算反应动力学速率常数的估计值
Figure DEST_PATH_IMAGE018
，得到无偏估计结果为
Figure DEST_PATH_IMAGE019
；
（3.3）重复步骤（3.2），将结果取平均，得到最终反应动力学速率常数无偏估计结果。
2.根据权利要求1所述反应动力学速率常数的无偏估计方法，其特征在于，步骤一具体为：
所述光谱数据矩阵D为：
Figure DEST_PATH_IMAGE020
，其中，C表示具有光吸收特性的物质随时间变化的浓度矩阵，S是每种物质在不同测量波长下的吸光系数矩阵，D表示反应溶液对于每种测量波长的吸光程度随时间的变化，
Figure 457841DEST_PATH_IMAGE003
与D维度相同；
以时间为自变量的常微分方程组建立反应动力学模型，给定光谱数据矩阵D和所有反应物的初始浓度，根据反应动力学速率常数的初始估计值
Figure 747002DEST_PATH_IMAGE011
，计算出C、S、
Figure 210345DEST_PATH_IMAGE003
，其中，C +为C的伪逆：
Figure DEST_PATH_IMAGE021
。
CN201910802085.8A 2019-08-28 2019-08-28 一种基于光谱的反应动力学速率常数的无偏估计方法 Active CN110514608B (zh)

Priority Applications (1)

Application Number Priority Date Filing Date Title
CN201910802085.8A CN110514608B (zh) 2019-08-28 2019-08-28 一种基于光谱的反应动力学速率常数的无偏估计方法

Applications Claiming Priority (1)

Application Number Priority Date Filing Date Title
CN201910802085.8A CN110514608B (zh) 2019-08-28 2019-08-28 一种基于光谱的反应动力学速率常数的无偏估计方法

Publications (2)

Publication Number Publication Date
CN110514608A CN110514608A (zh) 2019-11-29
CN110514608B true CN110514608B (zh) 2021-08-24

Family

ID=68628676

Family Applications (1)

Application Number Title Priority Date Filing Date
CN201910802085.8A Active CN110514608B (zh) 2019-08-28 2019-08-28 一种基于光谱的反应动力学速率常数的无偏估计方法

Country Status (1)

Country Link
CN (1) CN110514608B (zh)

Families Citing this family (1)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN112818558A (zh) * 2021-02-22 2021-05-18 浙江工业大学 一种针对反应过程机理建模的可估参数子集选择方法

Citations (9)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
US7181323B1 (en) * 2004-10-25 2007-02-20 Lockheed Martin Corporation Computerized method for generating low-bias estimates of position of a vehicle from sensor data
CN1916605A (zh) * 2006-08-31 2007-02-21 山东省科学院海洋仪器仪表研究所 化学发光动力学曲线分析臭氧氧化反应动力学参数的方法
CN101001219A (zh) * 2006-01-09 2007-07-18 电子科技大学中山学院 收敛速度快的特征参数自适应盲估计方法
CN101218242A (zh) * 2005-04-01 2008-07-09 孟山都技术公司 在草甘膦的生产中对n-(膦酰甲基)亚氨基二乙酸转化率的控制
CN101260037A (zh) * 2008-03-17 2008-09-10 华东理工大学 一种烷基芳烃液相氧化动力学机理模型的建模方法
CN101508768A (zh) * 2009-03-20 2009-08-19 华东理工大学 一种工业聚酯生产过程的智能建模方法
CN102744379A (zh) * 2012-03-07 2012-10-24 中冶南方工程技术有限公司 基于卡尔曼滤波的结晶器控制***状态估计方法
CN103105778A (zh) * 2013-02-05 2013-05-15 中原工学院 一种工业过程仿真数学模型参数的估计方法
CN109270484A (zh) * 2018-07-24 2019-01-25 南京航空航天大学 一种基于运动合成阵列的多信源doa估计方法

Family Cites Families (1)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
WO2015066052A1 (en) * 2013-10-28 2015-05-07 New York University Methods, computer-accessible medium and systems to model disease progression using biomedical data from multiple patients

Patent Citations (9)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
US7181323B1 (en) * 2004-10-25 2007-02-20 Lockheed Martin Corporation Computerized method for generating low-bias estimates of position of a vehicle from sensor data
CN101218242A (zh) * 2005-04-01 2008-07-09 孟山都技术公司 在草甘膦的生产中对n-(膦酰甲基)亚氨基二乙酸转化率的控制
CN101001219A (zh) * 2006-01-09 2007-07-18 电子科技大学中山学院 收敛速度快的特征参数自适应盲估计方法
CN1916605A (zh) * 2006-08-31 2007-02-21 山东省科学院海洋仪器仪表研究所 化学发光动力学曲线分析臭氧氧化反应动力学参数的方法
CN101260037A (zh) * 2008-03-17 2008-09-10 华东理工大学 一种烷基芳烃液相氧化动力学机理模型的建模方法
CN101508768A (zh) * 2009-03-20 2009-08-19 华东理工大学 一种工业聚酯生产过程的智能建模方法
CN102744379A (zh) * 2012-03-07 2012-10-24 中冶南方工程技术有限公司 基于卡尔曼滤波的结晶器控制***状态估计方法
CN103105778A (zh) * 2013-02-05 2013-05-15 中原工学院 一种工业过程仿真数学模型参数的估计方法
CN109270484A (zh) * 2018-07-24 2019-01-25 南京航空航天大学 一种基于运动合成阵列的多信源doa估计方法

Non-Patent Citations (3)

* Cited by examiner, † Cited by third party
Title
New Approach for Investigating Reaction Dynamics and Rates with Ab Initio Calculations;Kelly L. Fleming et al.;《The Journal of Physical Chemistry A》;20151221;第120卷;第299-305页 *
基于连串反应在线光谱监测数据的动力学曲线拟合及速率常数的估计;陈迪钊 等;《计算机与应用化学》;20020728;第19卷(第4期);第399-408页 *
气体No／s02／N2／02等离子体反应机理及动力学分析;李婷婷 等;《工程热物理学报》;20070531;第28卷(第3期);第532-533页 *

Also Published As

Publication number Publication date
CN110514608A (zh) 2019-11-29

Similar Documents

Publication Publication Date Title
Benhamou et al. Time dependent Heston model
Kelleher et al. Atomic transition probabilities of sodium and magnesium. A critical compilation
Wu et al. The QCD renormalization group equation and the elimination of fixed-order scheme-and-scale ambiguities using the principle of maximum conformality
Noor et al. Consistent estimation of Gibbs energy using component contributions
Bozzi et al. Production of Drell–Yan lepton pairs in hadron collisions: Transverse-momentum resummation at next-to-next-to-leading logarithmic accuracy
Royston et al. Approximating statistical functions by using fractional polynomial regression
Bhatt et al. Incremental identification of reaction systems—A comparison between rate-based and extent-based approaches
Zaika et al. Parametric identification of hydrogen permeability model by delay times and conjugate equations
Yurchenko et al. Ab initio dipole moment and theoretical rovibrational intensities in the electronic ground state of PH3
CN110514608B (zh) 一种基于光谱的反应动力学速率常数的无偏估计方法
Kowzan et al. Fully quantum calculations of the line-shape parameters for the Hartmann-Tran profile: A CO-Ar case study
Jaumot et al. Noise propagation and error estimations in multivariate curve resolution alternating least squares using resampling methods
Colli et al. Well-posedness and optimal control for a Cahn-Hilliard-Oono system with control in the mass term
Chen et al. Kinetic parameter estimation based on spectroscopic data with unknown absorbing species
Lazopoulos et al. NLO QCD corrections to the production of tt¯ Z in gluon fusion
Gracey Four loop MS mass anomalous dimension in the Gross–Neveu model
Feng et al. Geometric Brownian motion with affine drift and its time-integral
Lesnoff et al. Monte Carlo methods for estimating Mallows's Cp and AIC criteria for PLSR models. Illustration on agronomic spectroscopic NIR data
Asaithambi Numerical solution of the one-dimensional time-independent Schrödinger’s equation by recursive evaluation of derivatives
Knoll et al. The thermodynamic scaling function of a polymer solution
Zogg et al. Identification of kinetic and thermodynamic reaction parameters from online calorimetric and IR-ATR data using a new combined evaluation algorithm
Naseri et al. A second-order standard addition method based on the data treatment by calculation of variation matrix of kinetic systems analyzed by MCR-ALS
Messick et al. Selecting factors for partial least squares
CN114334030A (zh) 基于量子支持向量机的高分子聚合反应产物评价方法
Bognarova et al. A comparative study of two numerical approaches for solving Kim–Kim–Suzuki phase-field models

Legal Events

Date Code Title Description
PB01 Publication
PB01 Publication
SE01 Entry into force of request for substantive examination
SE01 Entry into force of request for substantive examination
GR01 Patent grant
GR01 Patent grant