2021年国赛高教杯数学建模B题乙醇偶合制备C4烯烃解题全过程文档及程序
2021年国赛高教杯数学建模
B题 乙醇偶合制备C4烯烃
原题再现
C4 烯烃广泛应用于化工产品及医药的生产，乙醇是生产制备 C4 烯烃的原料。在制备过程中，催化剂组合（即：Co 负载量、Co/SiO2 和 HAP 装料比、乙醇浓度的组合）与温度对 C4 烯烃的选择性和 C4 烯烃收率将产生影响（名词解释见附录)。因此通过对催化剂组合设计，探索乙醇催化偶合制备 C4 烯烃的工艺条件具有非常重要的意义和价值。
某化工实验室针对不同催化剂在不同温度下做了一系列实验，结果如附件 1 和附件 2 所示。请通过数学建模完成下列问题：
(1) 对附件 1 中每种催化剂组合，分别研究乙醇转化率、C4 烯烃的选择性与温度的关系，并对附件 2 中 350 度时给定的催化剂组合在一次实验不同时间的测试结果进行分析。
(2) 探讨不同催化剂组合及温度对乙醇转化率以及 C4 烯烃选择性大小的影响。
(3) 如何选择催化剂组合与温度，使得在相同实验条件下 C4 烯烃收率尽可能高。若使温度低于 350 度，又如何选择催化剂组合与温度，使得 C4 烯烃收率尽可能高。
(4) 如果允许再增加 5 次实验，应如何设计，并给出详细理由。
整体求解过程概述(摘要)
本文主要研究乙醇偶合制备 C4 烯烃时，不同催化剂组合及温度对乙醇转化率和 C4 烯经选择性的影响、C4 烯经收率最高时的催化剂组合与温度最优化的实验设计问题。
针对问题一，分析附件 1 中每种催化剂组合数据规律,分别建立乙醇转化率、C4 烯烃选择性与温度之间的一元二次回归模型。基于最小二乘法,运用 Matlab 中的 polyfit()和 regress()函数对该模型进行求解，得到回归系数，然后通过判别系数 R^2 判断模型是否准确，如果判别系数小于 0.6，则认为该模型不准确，需要进行进一步通过剔除离总体偏差较大的数据与优化方程模型进行模型的优化，进而通过分析所得曲线的变化规律得到结论:在一定范围内随着温度升高，乙醇的转化率与 C4 烯怪的选择性均会增高，但若温度过高，则 C4 烯经的选择性可能会降低；对附件 2 中的数据进行对比分析，并结合化学反应机理得到结论:乙醇转化率随着时间增大逐渐减小，C4 烯烃选择性随时间变化的趋势较为平稳，随着反应时间增加，各产物的选择性都趋于平稳，即化学反应达到平衡状态。
针对问题二，需要分析不同催化剂组合及温度对乙醇转化率和 C4 烯烃选择性大小的影响。首先，用控制变量的方法对附件 1 中的数据进行分类对比，通过定性分析得到相应结论。得到: Co 负载量对 C4 烯烃选择性的影响最大，Co/SiO2 和 HAP 装料比对 C4 烯烃选择性的影响最小。最后，结合实际考虑因素间交互作用，分别建立 C4 烯烃选择性、乙醇转化率与各因素及其交互作用的多元线性回归模型，通过对问题二进行比较分析法得出：Co 负载量与Co/SiO2 和 HAP 装料比共同作用对乙醇转化率影响最大，温度对 C4 烯烃的选择性影响最大。
针对问题三，首先在问题二的基础上进行数据预处理，除去 A11 组有石英砂无 HAP 的一组数据并在问题二多元回归模型的基础上，给出 C4 烯烃收率的函数表达式。其次以温度、Co/SiO2 和 HAP 的质量和、乙醇浓度、Co 负载量作为决策变量，以最高 C4 烯烃收率作为目标函数，以各个影响因素的取值范围作为约束条件,建立 C4 烯烃收率最高的优化模型。然后基于迭代算法，利用最小二乘法拟合，然后运用牛顿迭代法对该模型求解，然后用选择排序法得到: C4 烯烃最优收率为 53%。再改变温度约束条件，使得最高温度不超过 350 度，得到此时 C4 烯烃最优收率为 12%。
针对问题四，结合前三问的分析结论，通过过前面问题的研究不难得出当温度达到一定时候，而对前三问的分析可以明显得出当温度为一范围内，温度越高 C4 烯烃转化率反而越低，C4 烯烃收率将会显著下降到一个最低点，而过了这一段范围后，C4 烯烃转率又随着温度增高而增高所以这里建立方程模型得出这个最低点与其变化规律，所以通过其他条件不变，温度变化，得出 C4 烯烃收率是否与温度的是否成二次曲线模型变化的二次函数方程模型。
模型假设：
1. 忽略外界因素对实验的影响
2. 忽略有石英砂无 HAP 对实验的影响
3. 环境温度和压强对实验的进行没有影响
4. 实验设备完好不会导致实验结果出错或实验失败
5. 判别系数 R^2>0.6 认为拟合有效
问题分析：
问题一的分析
问题一分别分析温度对乙醇转化率、C4 烯烃选择性的关系。问题一的第一部分，首先对附件 1 画出每种催化剂组合的乙醇转化率、C4 烯烃选择性与温度的散点图，在此基础上我们初步得出乙醇转化率总体上随着温度的升高而增大的结论，但是某些数值在增加到某个温度时会有所下降，我们需要对异常的点进行处理，然后进行分类。如果乙醇转化率、C4 烯烃选择性与温度具有一定的线性关系，则建立数学模型利用MATLAB 求出相关系数。
问题一的第二部分，我们根据图像拟合对化学反应的进程和产物进行分析，得出在 350°时乙醇转化率、C4 烯烃选择性、C4 烯烃收率等指标随时间的变化规律，分析出反应过程中的产物量比及反应的平衡时间。
问题二的分析
问题二是分析不同催化剂组合及温度对乙醇转化率和 C4 烯烃选择性大小的影响。首先我们通过控制变量的方法，对 20 组试验数据进行组合对照，进行直观的定量分析。控制自变量 Co 负载量（Co/SiO2 的重量之比）、Co/SiO2 和 HAP 装料比Co/SiO2 和 HAP 质量之和，乙醇浓度、温度对因变量乙醇转化率、C4 烯烃选择性的影响。其间我们需要发现自变量不同催化剂组合、温度之间的重要因素和次要因素，并在定性分析的基础上进行定量分析，确定变量之间的关系，找出合适的数学表达式。
问题三的分析
基于问题二的模型建立，首先对数据进行预处理，去除 A11 组有石英砂无 HAP 这一组的实验数据，对 C4 烯烃收率的函数表达式以温度、Co/SiO2 和 HAP 的质量和、乙醇浓度、Co负载量，Co/SiO2 和 HAP 装料比作为决策变量，以最高 C4 烯烃收率作为目标函数，以各个影响因素的取值范围作为约束条件,建立 C4 烯烃收率最高的优化模型。
问题四的分析
综合考虑问题一、二、三,对于附件中没有改变某个实验条件可能使得实验结果更优的进行实验，通过对前三问的分析可以明显得出当温度为一范围内，温度越高 C4 烯烃转化率反而越低，C4 烯烃收率将会显著下降到一个最低点，而过了这一段范围后，C4 烯烃转率又随着温度增高而增高。但是对于因缺少数据而使得某个结论缺服说服力的进行实验。
模型的建立与求解整体论文缩略图
全部论文请见下方“ 只会建模 QQ名片” 点击QQ名片即可
程序代码：(代码和文档not free)
%对数据进行导入
%导入附件 1 中的数据
%wendu - 温度 yczhl - 乙醇转化率 wyl - 对温度与乙醇转化率进行拟合后得到的函数
wenduA1=xlsread("C:UsersASUSDesktop 数 学 建 模 题 目 B 附 件 1.xlsx"," 性 能 数 据 表
","C2:C6");
yczhlA1=xlsread("C:UsersASUSDesktop数学建模题目B附件 1.xlsx","性能数据表","D2:D6");
wenduA2=xlsread("C:UsersASUSDesktop 数 学 建 模 题 目 B 附 件 1.xlsx"," 性 能 数 据 表
","C7:C11");
yczhlA2=xlsread("C:UsersASUSDesktop 数 学 建 模 题 目 B 附 件 1.xlsx"," 性能数据表
","D7:D11");
wenduA3=xlsread("C:UsersASUSDesktop 数 学 建 模 题 目 B 附 件 1.xlsx"," 性 能 数 据 表
","C12:C18");
yczhlA3=xlsread("C:UsersASUSDesktop 数 学 建 模 题 目 B 附 件 1.xlsx"," 性能数据表
","D12:D18");
wenduA4=xlsread("C:UsersASUSDesktop 数 学 建 模 题 目 B 附 件 1.xlsx"," 性 能 数 据 表
","C19:C24");
yczhlA4=xlsread("C:UsersASUSDesktop 数 学 建 模 题 目 B 附 件 1.xlsx"," 性能数据表
","D19:D24");
wenduA5=xlsread("C:UsersASUSDesktop 数 学 建 模 题 目 B 附 件 1.xlsx"," 性 能 数 据 表
","C25:C30");
yczhlA5=xlsread("C:UsersASUSDesktop 数 学 建 模 题 目 B 附 件 1.xlsx"," 性能数据表
","D25:D30");
wenduA6=xlsread("C:UsersASUSDesktop 数 学 建 模 题 目 B 附 件 1.xlsx"," 性 能 数 据 表
","C31:C35");
%k 为拟合次数
k=3;
wyA1=polyfit(wenduA1,yczhlA1,k); %求出拟合参数
A1=vpa(poly2sym(wyA1),10);%求出曲线 %画出拟合曲线为
y1=polyval(wyA1,wenduA1);
plot(wenduA1,y1,wenduA1,yczhlA1,"*r");
hold on ;
wyA2=polyfit(wenduA2,yczhlA2,k);
A2=vpa(poly2sym(wyA2),10);%曲线为
%画出拟合曲线为
y2=polyval(wyA2,wenduA2);
plot(wenduA2,y2,wenduA2,yczhlA2,"*r");
hold on
wyA3=polyfit(wenduA3,yczhlA3,k);
A3=vpa(poly2sym(wyA3),10);%曲线为
%画出拟合曲线为
y3=polyval(wyA3,wenduA3);
plot(wenduA3,y3,wenduA3,yczhlA3,"*r");
hold on ;
wyA4=polyfit(wenduA4,yczhlA4,k);
A4=vpa(poly2sym(wyA4),10);%曲线为
%画出拟合曲线为
y4=polyval(wyA4,wenduA4);
plot(wenduA4,y4,wenduA4,yczhlA4,"*r");
hold on
wyA5=polyfit(wenduA5,yczhlA5,k);
A5=vpa(poly2sym(wyA5),10);%曲线为
%画出拟合曲线为
y5=polyval(wyA5,wenduA5);
plot(wenduA5,y5,wenduA4,yczhlA5,"*r");
hold on
wyA6=polyfit(wenduA6,yczhlA6,k);
A6=vpa(poly2sym(wyA6),10);%曲线为
%画出拟合曲线为
y6=polyval(wyA6,wenduA6);
plot(wenduA6,y6,wenduA6,yczhlA6,"*r");
hold on
wyA7=polyfit(wenduA7,yczhlA7,k);
A7=vpa(poly2sym(wyA7),10);%曲线为
%画出拟合曲线为
y7=polyval(wyA7,wenduA7);
plot(wenduA7,y7,wenduA7,yczhlA7,"*r");
hold on
wyB7=polyfit(wenduB7,yczhlB7,k);
B7=vpa(poly2sym(wyB7),10);%曲线为
%画出拟合曲线为
y21=polyval(wyB7,wenduB7);
plot(wenduB7,y21,wenduB7,yczhlB7,"*r");
%{
legend("A1=0.3331893326*x - 84.08276818","A2=0.6629580734*x -
161.8909145","A3=0.4195648852*x - 95.88308918",...
"A4=0.5817118467*x - 144.5709001","A5=0.4078621995*x -
97.62301858","A6=0.5015334253*x - 119.8325796",...
"A7=0.3775367296*x - 74.26046445","A8=0.3396341323*x - 83.77607016",
"A9=0.2491040369*x - 65.67280331",...
"A10=0.1834135666*x - 49.64860852","A11=0.2067687555*x -
56.50163086","A12=0.2859331347*x - 74.81353909",...
"A13=0.3359113498*x - 86.68673211","A14=0.3359113498*x -
86.68673211","B1=0.2794712209*x - 73.14252889",...
"B2=0.2726886201*x - 70.9546115","B3=0.1313117907*x -
36.1563065","B4=0.2077582601*x - 56.81340574","B5=0.2715966155*x - 72.27286731",...
"B6=0.3837143453*x - 100.0040718","B7=0.4197142857*x - 109.3428571");
%}
legend('A1','A2','A3','A4','A5','A6','A7','A8','A9','A10','A11','A12','A14','A14','B1','B2','B3','B4','B5','B
6','B7')
xlabel("温度");
ylabel("乙醇转化率");
title("温度-乙醇转化率拟合图")
avgA1=(sum(yczhlA1))/length(yczhlA1);%乙醇转化率平均值
totA1=sum((yczhlA1-avgA1).^2);
resA1= double( sum(( subs(A1,wenduA1)-yczhlA1).^2));
R2A1=1-resA1/totA1;
avgA2=(sum(yczhlA2))/length(yczhlA2);%乙醇转化率平均值
totA2=sum((yczhlA2-avgA2).^2);
resA2= double(sum((subs(A2,wenduA2)-yczhlA2).^2));
R2A2=1-resA2/totA2;