gPROMS二次开发教程（15）：优化与最优控制——动态优化、求解器选项与结果解读

gPROMS二次开发教程15优化与最优控制——动态优化、求解器选项与结果解读版本声明块工具/软件gPROMS 桌面建模环境gPROMS ModelBuilder检索期官方发布锚点gPROMS Process 2022.1.0适用版本以官方发布说明为准语言/环境gPROMS 语言equation-oriented modelling结果后处理脚本用 Python 3 NumPy离线复算不调用任何 gPROMS 组件本文目标把模型能算推进到模型能决策——学会定义优化问题、读懂优化结果并且不信优化结果一句话结论优化就是在约束围成的可行域里找目标的最小/最大值官方能力线为Steady state flowsheeting and optimization与Dynamic flowsheeting and optimization参数估计Estimating model parameters from experimental data本质上是同一套优化引擎的另一种用法而优化结果在通过约束活性 多初值复跑 可行性复核三项体检之前一律不能进报告。〇、本篇要解决的认知问题Q1gPROMS 的优化能力定位是什么为什么说参数估计和优化是同一台机器Q2优化问题必须写清哪三要素少写一个会发生什么Q3动态优化与稳态优化在建模与结果形态上有什么本质差别Q4求解器选项与收敛策略该怎么设为什么本文只给边界不给经验值Q5优化结果为什么要做可信度检查局部最优、约束活性、初值依赖分别怎么查一、机制解析1.1 gPROMS 的优化能力定位gPROMS 以equation-oriented modelling方程导向建模为核心特征模型是开放式白箱方程组因此优化不必像黑箱流程模拟器那样靠外挂迭代而是可以把决策变量直接交给平台、在方程组层面统一求解与搜索。官方资源目录中的优化相关能力名包括Steady state flowsheeting and optimization、Dynamic flowsheeting and optimization与Estimating model parameters from experimental data。关键认知参数估计就是一次特殊的优化——目标函数是加权残差平方和约束是参数取值范围因此初值敏感、解不唯一、不可辨识会在优化语境里原样复现为初值依赖与局部最优而优化器不创造物理它只在你给的模型与约束内寻优。1.2 优化三要素决策变量 / 目标 / 约束要素定义少写的后果决策变量允许被搜索的量进料温度、冷却水流量、温控曲线节点值、参数值会去优化不该动的量如设备尺寸得到无法实施的方案目标函数单一标量能耗最小、收率最大、时间最短、残差最小多目标混写“又省能耗又好质量”会让求解器无法比较优劣约束等式约束模型方程不等式约束安全/质量/设备上限丢掉不等式约束 → 得到数学最优、物理危险的解如超温┌──────────────┬────────────────────────────────────────────┐ │ 决策变量 │ u [T_feed, F_cw, 温控曲线节点 T_c(t_i)] │ │ 目标函数 │ min J 冷却能耗 E_cw或残差平方和 SSE │ │ 等式约束 │ 模型方程物料/能量衡算、速率式 │ │ 不等式约束 │ T_reactor(t) ≤ 400 KX ≥ 0.90F_cw ≥ 0 │ │ 自由度核验 │ 第 1 次模型自身自由度 0能解 │ │ │ 第 2 次放开的决策变量数 可搜索自由度 │ └──────────────┴────────────────────────────────────────────┘铁律 8建模先算自由度优化前必须做两次数自由度——第一次要求模型自身自由度为 0模型能解第二次要求放开的量恰好等于决策变量个数优化有东西可搜。两次都对不上优化一定失败。1.3 动态优化 vs 稳态优化两条官方能力线的差别不在算法而在问题形态维度Steady state flowsheeting and optimizationDynamic flowsheeting and optimization决策变量一组常数操作条件、设备参数一条轨迹或其若干节点值约束代数不等式如T ≤ 400 K路径约束全程满足如T(t) ≤ 400 K ∀t与终点约束结果形态一个最优点一条最优曲线 对应状态轨迹典型用途稳态流程能耗/收率优化、设计裕度最优温控曲线、开车/停车策略、批次配方优化动态优化的三条硬约束①路径约束是全程约束只在节点上检查会漏掉节点之间的超温②决策变量的节点数就是曲线分辨率它同时决定问题规模与曲线能否平滑可实施③最优曲线必须能落到现场执行——若最优解要求每分钟翻转一次阀门那不是最优控制那是数学游戏。经验法则动态优化的输出永远先过可实施性审查变化率是否落在执行机构能力内、曲线是否平滑。不可实施的最优 不可用。1.4 求解器选项与收敛策略为什么只给边界优化求解器涉及算法选择、收敛判据、可行性容差、迭代上限等大量选项它们属于闭源实现细节✅ 可以说存在求解器选项与收敛相关设置调参前先做问题体检不可行 / 病态 / 多解 ❌ 不可以说某选项的默认值是多少、某容差该设成 1e-6 还是 1e-8、某算法一定更快本篇给出的是与闭源细节无关的通行纪律在任何优化器上都成立纪律做法理由先可行、再最优先找一个满足全部约束的解从不可行起点出发优化器可能绕远路或直接失败无量纲化把量级差异大的变量/约束做尺度归一条件数差的优化问题收敛极慢、结果不可信先简化、后加密先用少量决策变量跑通再加密曲线节点便于区分模型问题与规模问题宁可多跑几组初值至少 3 组不同初值单次收敛不能证明是全局最优见 1.5 节铁律 4商业闭源软件的诚实边界求解器选项清单、默认值、算法名称与调优建议一律以官方文档为准本文不做任何数值层面的臆造。1.5 优化结果的可信度检查三项体检体检项现象怎么查处置局部最优换个初值结果明显更好多初值复跑≥3 组比目标值报告多个解的差异必要时改初值策略约束活性最优解压在约束边界上检查哪些不等式约束取等号间隙为 0明确写出受约束控制若无意顶边说明约束写漏或过紧初值依赖目标值随初值漂移记录初值 → 目标对照表把初值写进结果记录敏感的解要标注风险另有两个必须养成的动作可行性复核把最优解回代原模型逐条核对不等式约束不采信求解器自报的可行与物理合理性复核负流量、超压、超温一律否决。铁律 9改完必验优化报告必须包含初值 → 目标值对照表与约束活性清单只给一条最优曲线和一句优化完成的报告等同于没做优化。1.6 把优化嵌进工作流优化天然是批量任务多初值、多工况初值生成用外部编排脚本gO:Python侧循环、逐场景求解并落盘为固定字段 CSV便于与回归基线比对接第 13 篇、必要时用gO:MATLAB与控制器算法闭环。二、完整代码/模型与逐行剖析代码 2-1非等温 CSTR 动态模型 温控曲线参数化结构示意% % 模型名 : cstr_noniso_dyn 模型类型 : LUMPED集中参数 / 动态 % 物理背景 : 非等温 CSTRA - B一级不可逆、放热用冷却夹套移热 % 温控曲线由若干节点值参数化供动态优化搜索。 % 说明 : 结构示意。PARAMETER / VARIABLE / LUMPED 为已确证关键字时间导数的 % 书写方式与把某量声明为决策变量的语法均以官方文档为准。 % 数据来源 : 数值为教学示意值正式使用须替换为可引用文献值铁律 10 % MODELS PARAMETER R 8.314 % 通用气体常数 [J/(mol*K)]文献值 V 1.0e-3 % 反应器有效体积 [m^3] F 1.0e-4 % 进料体积流量 [m^3/s] C_Ain 1.0 % 进料浓度 [kmol/m^3] T_in 300.0 % 进料温度 [K] k0 8.0e6 % 指前因子 [1/s] Ea 5.8e4 % 活化能 [J/mol] rho_cp 2.0e6 % 体积热容 [J/(m^3*K)] dH -8.0e7 % 反应热 [J/kmol]放热取负 UA 2.0e3 % 总传热系数×面积 [J/(s*K)] T_c1 290.0 % 冷却介质温度节点 1 [K] ← 决策变量示意 T_c2 290.0 % 冷却介质温度节点 2 [K] ← 决策变量示意 T_c3 290.0 % 冷却介质温度节点 3 [K] ← 决策变量示意 VARIABLE C_A % A 的浓度 [kmol/m^3]状态量 T % 反应器温度 [K]状态量 r % 反应速率 [kmol/(m^3*s)]辅助量 k % 速率常数 [1/s]辅助量 T_c % 冷却介质温度 [K]由节点插值得到 dCA_dt % C_A 对时间的一阶导数项记号以官方文档为准 dT_dt % T 对时间的一阶导数项记号以官方文档为准 %% 方程区残差式 k - k0*exp(-Ea/(R*T)) 0 % 本构Arrhenius r - k*C_A 0 % 本构一级反应速率 dCA_dt r - F*(C_Ain - C_A)/V 0 % 动态物料衡算残差 rho_cp*dT_dt - r*(-dH) - F*rho_cp*(T_in - T)/V UA*(T - T_c)/V 0 %% ↑ 蓄热项 ↑ 反应放热 ↑ 进料带入/带出 ↑ 冷却夹套移热 %% 量纲核对各项均为 [J/(m^3*s)]等式两侧同量纲 %% 时间导数dCA_dt、dT_dt与节点插值关系的写法以官方文档为准逐行剖析T_c1~T_c3放在PARAMETER段但在优化中被放开与第 11 篇待估参数也是参数同构——被优化的量依然是参数只是交给了优化流程去搜索1.2 节。把温控曲线参数化成有限个节点是动态优化最关键的建模动作它把搜索一条函数降级为搜索几个数。两条衡算都是残差式冷却项UA*(T - T_c)/V显式出现这一点在优化中极其重要——若目标函数含冷却能耗它与模型的耦合点就是这一项。能量衡算里看不到T_c优化器就根本动不到这个量。自由度核验铁律 8待求量 7 个C_A、T、r、k、T_c、dCA_dt、dT_dt本块给出 4 条方程余下约束来自时间导数定义与插值关系。账必须算平每条待求量都要有对应方程否则模型不可解、优化无从谈起。代码 2-2优化问题定义表决策变量 / 目标 / 约束% % 优化问题定义表动态优化非等温 CSTR 最优温控曲线—— 本表是优化问题的规格书 % 【问题类型】 Dynamic flowsheeting and optimization 官方能力名 对照决策变量为常数、约束为代数不等式时用 Steady state flowsheeting and optimization 【决策变量】 u1..u3 T_c1..T_c3 冷却介质温度节点 [K] 范围 [280, 320] 节点数即曲线分辨率先少后多见 1.4 节 【目标函数】 min J ∫ |UA*(T - T_c)| dt % 与模型冷却项直接对应 【约束】 等式Arrhenius 速率式、一级速率式、动态物料/能量衡算 不等式路径约束全程满足T(t) ≤ 400 KT(t) ≥ 285 K X_final ≥ 0.90T_c_i ∈ [280, 320] K 【实施性约束】dT_c/dt 变化率落在执行机构能力内相邻节点差不宜过大曲线平滑 【自由度核验】第 1 次模型自身自由度 0第 2 次决策变量数 3铁律 8 【结果体检】 ① 多初值复跑比对目标值 ② 列出取等号的活性约束 ③ 最优解回代复核不等式约束 ④ 物理合理性复核逐行剖析目标函数与模型的耦合点被点出J ∫ |UA*(T - T_c)| dt对应代码 2-1 的冷却项这是优化最常见的失败原因之一——目标里的量与模型里的决策变量不是同一个优化器会优化了个寂寞目标值几乎不变。不等式约束显式标注为路径约束落地 1.3 节硬约束①实施性约束单独成栏把最优曲线能否落地执行写进问题定义而不是留给事后争论。代码 2-3稳态流程能耗最小化结构示意% % 模型名 : cstr_ss_energy_opt 能力对应 : Steady state flowsheeting and optimization % 物理背景 : 稳态 CSTRA - B放热用冷却水移热目标为最小化冷却能耗。 % 决策变量 : T_feed进料温度、F_cw冷却水流量 % 说明 : 结构示意段名与关键字以官方文档为准数值为教学示意值铁律 10 % 自由度 : 待求量 7 个C_A_out, T, r, k, T_cw_out, X, E_cwvs 方程 7 条 ⇒ 0 % MODELS PARAMETER R 8.314 % 通用气体常数 [J/(mol*K)]文献值 V 1.0e-3 % 反应器有效体积 [m^3] F 1.0e-4 % 进料体积流量 [m^3/s] C_Ain 1.0 % 进料浓度 [kmol/m^3] k0 8.0e6 % 指前因子 [1/s] Ea 5.8e4 % 活化能 [J/mol] rho_cp 2.0e6 % 体积热容 [J/(m^3*K)] dH -8.0e7 % 反应热 [J/kmol] cp_cw 4.18e6 % 冷却水体积热容 [J/(m^3*K)] T_feed 300.0 % 进料温度 [K] ← 决策变量 1范围 [295, 320] F_cw 5.0e-6 % 冷却水流量 [m^3/s] ← 决策变量 2范围 [1e-6, 5e-5] VARIABLE C_A_out % 出口浓度 [kmol/m^3] T % 反应器温度 [K] r % 反应速率 [kmol/(m^3*s)] k % 速率常数 [1/s] T_cw_out % 冷却水出口温度 [K] X % 转化率 [—]受约束的关键输出 E_cw % 冷却能耗 [J/s]目标函数示意定义 %% 方程区残差式四条主方程 三条定义式 k - k0*exp(-Ea/(R*T)) 0 % 本构Arrhenius r - k*C_A_out 0 % 本构一级速率 F*(C_Ain - C_A_out) - V*r 0 % 守恒物料衡算 F*rho_cp*(T_feed - T) - V*r*(-dH) - F_cw*cp_cw*(T_cw_out - T) 0 %% ↑ 进料带入/带出 ↑ 反应放热 ↑ 冷却水移热含决策变量 F_cw X - (C_Ain - C_A_out)/C_Ain 0 % 定义转化率 T_cw_out - (T - 5.0) 0 % 简化假设冷却水温差 5 K E_cw - F_cw*cp_cw*(T_cw_out - 298.15) 0 % 定义冷却能耗逐行剖析两个决策变量都在能量衡算里出现T_feed在F*rho_cp*(T_feed - T)中、F_cw在F_cw*cp_cw*(T_cw_out - T)中。决策变量必须在模型方程里看得见——否则它不影响任何结果搜索它只是浪费算力1.2 节耦合点纪律。X与E_cw显式声明为变量X是受约束的输出X ≥ 0.90E_cw是目标函数。把目标与约束都变成可寻址的名字优化配置才有抓手。T_cw_out - (T - 5.0) 0是显式标注的简化假设把它写成方程而不是心算掉将来换成真实换热器模型时只需替换这一条其余不动——假设也要显式化是长模型可维护性的核心技巧。代码 2-4优化结果体检脚本约束活性 多初值复跑 可行性复核# -*- coding: utf-8 -*- 优化结果体检脚本对优化器给出的最优解做三项独立检查。 纪律不采信求解器自报的可行/最优一律回代原模型逐条复核 多初值复跑的初值必须写进结果记录本脚本为离线复算不调用 gPROMS 组件。 importnumpyasnp# --- 与代码 2-3 一致的模型参数复核必须回到原模型 -------------------------R_GAS,V,F,C_AIN8.314,1.0e-3,1.0e-4,1.0K0,EA,RHO_CP,DH8.0e6,5.8e4,2.0e6,-8.0e7CP_CW4.18e6defsteady_state(T_feed,F_cw):给定两个决策变量解出稳态输出此处用解析近似式复核 真实场景应把最优解回代 gPROMS 模型求解一次取回同一组输出。TT_feed40.0*(F_cw/5.0e-6)**(-0.35)# 冷却越强温度越低示意kK0*np.exp(-EA/(R_GAS*T))# Arrhenius 速率常数C_A_outF*C_AIN/(FV*k)# 稳态出口浓度解析解X(C_AIN-C_A_out)/C_AIN# 转化率E_cwF_cw*CP_CW*5.0# 冷却能耗示意定义return{T:T,X:X,C_A_out:C_A_out,E_cw:E_cw}defcheck_feasibility(T_feed,F_cw,X_min0.90,T_max400.0,T_min285.0):可行性复核逐条核对不等式约束返回明细。不采信求解器自报结论。ssteady_state(T_feed,F_cw)checks[(T 400 K,s[T]T_max,fT {s[T]:.2f}K),(T 285 K,s[T]T_min,fT {s[T]:.2f}K),(X 0.90,s[X]X_min,fX {s[X]:.4f}),(F_cw 0,F_cw0.0,fF_cw {F_cw:.3e}m^3/s),]returnall(c[1]forcinchecks),checks,sdefcheck_activity(T_feed,F_cw,tol1.0e-3):约束活性检查判断哪些不等式约束取等号即被约束顶住。ssteady_state(T_feed,F_cw)gaps{T 400 K:400.0-s[T],X 0.90:s[X]-0.90}active[nameforname,gapingaps.items()ifabs(gap)tol]returngaps,activeif__name____main__:# ① 多初值复跑汇总查局部最优初值由第 14 篇的编排骨架批量驱动runs[(初值A 295K/1e-6,295.0,1.0e-6),(初值B 305K/1e-5,305.0,1.0e-5),(初值C 315K/3e-5,315.0,3.0e-5)]print( ① 多初值复跑汇总查局部最优)fortag,tf,fcinruns:ok,_,scheck_feasibility(tf,fc)print(f{tag}: E_cw {s[E_cw]:.4e}J/s 可行{ok})print( 判读若三组目标值差异显著说明存在局部最优须择优并标注风险。)# ② 约束活性 ③ 可行性复核只对当前最优解逐条核对best(T_feed 305.0 K, F_cw 1.0e-5,305.0,1.0e-5)print(\n ② 约束活性检查 )gaps,activecheck_activity(best[1],best[2])forname,gapingaps.items():print(f{name}: 间隙 {gap:.4f})print(f 活跃约束被顶住{activeifactiveelse无})print(\n ③ 可行性复核不采信求解器自报结论)ok,checks,scheck_feasibility(best[1],best[2])forname,passed,detailinchecks:print(f [{PASSifpassedelseFAIL}]{name:12s}{detail})print(整体可行性 ,PASSifokelseFAIL)逐行剖析check_feasibility返回名称, 是否通过, 实测值三元组实测值是灵魂——它让复核结论可审计“你说满足满足到什么程度”。只打印 PASS/FAIL 的复核等于没复核。check_activity输出间隙而非布尔值间隙为 0 说明该约束正顶住最优解。工程上极重要——若X ≥ 0.90恰好活跃说明转化率卡着下限任何扰动都会让它不合格若它完全不活跃则可能是被别的约束挡住了可考虑放宽那条再试。多初值复跑在脚本里是汇总而非重跑重跑由第 14 篇的编排骨架批量驱动本脚本负责汇总与判读——优化 批量编排正是第 14、15 两篇的接口。三、常见报错与排查报错 3-1优化收敛了但目标值几乎没变。现象优化前后结果一致。根因决策变量未出现在模型方程中耦合点断裂。解法在模型里搜索每个决策变量的名字确认它真的出现在某条方程里再确认目标函数含该变量。报错 3-2优化得到物理上不可行的解。现象出现负流量、超温、超压。根因不等式约束漏写或过宽三要素缺一件。解法补齐路径约束与实施性约束代码 2-2 的规格书。报错 3-3换个初值结果就好了很多。现象目标值随初值大幅漂移。根因局部最优。解法多初值复跑≥3 组记录初值 → 目标值对照表报告多个解的差异并标注风险。报错 3-4优化中途失败报不可行或残差不降。现象迭代无法继续。根因起始点不可行或量纲差异导致病态或模型自身自由度不为 0。解法先找可行起点“先可行、再最优”做无量纲化回头做第一次自由度核验。报错 3-5最优温控曲线太陡现场执行不了。现象数学最优但不可实施。根因未加实施性约束或节点过多导致曲线振荡。解法把dT_c/dt上界写进约束减少节点数并采用平滑参数化。铁律 3版本先声明求解器与优化器的选项名称、默认设置、可用算法在不同版本间存在差异以官方发布说明与随附文档为准。四、动手练习练习 1规格书落地用代码 2-2 的格式为自己的模型写一份完整优化问题规格书。判定决策变量、目标函数、等式约束、不等式约束、实施性约束、两次自由度核验、结果体检七栏齐全。练习 2耦合点检查把某个决策变量从能量衡算里删除仅用于实验重跑优化。判定观察到目标值不再改善能说明耦合点断裂为何导致优化失效随后改回。练习 3局部最优用三组不同初值做一次优化可借第 14 篇的编排骨架。判定给出初值 → 目标值对照表明确指出是否存在局部最优并给出择优依据。练习 4约束活性运行代码 2-4 的第二段。判定列出所有取等号的活跃约束及其间隙解释转化率下限恰好活跃对工艺裕度的含义。五、小结与下一篇预告本篇把模型从能算推到能决策官方两条能力线为Steady state flowsheeting and optimization与Dynamic flowsheeting and optimization而参数估计本质是同一台优化引擎优化问题必须写清三要素决策变量 / 目标 / 约束其中问题类型一栏的判据是决策变量为常数还是轨迹并做两次自由度核验动态优化的三条硬约束是路径约束全程有效、曲线分辨率即代价、最优曲线必须可实施求解器选项属闭源细节一律以官方文档为准本篇只给与实现无关的通行纪律结果必须过三项体检——多初值复跑、约束活性间隙、回代可行性复核。没有体检的优化结论等同于一条漂亮但不可信的曲线。第 16 篇《与其他工具与控制系统联用》把模型放进更大的系统CAPE-OPEN Unit Operation的向导式导出、SIMIT gPROMS coupling联合仿真原理 / 启用 / 配置、FMI/FMU作为跨工具载体以及与 CFD 的多区混合模型。本篇认知问题回显FAQQ1gPROMS 的优化能力定位是什么为什么说参数估计和优化是同一台机器AgPROMS 以方程导向建模为核心模型是开放式白箱方程组优化可直接在方程组层面统一求解与搜索。官方优化相关能力线包括稳态流程与优化Steady state flowsheeting and optimization、动态流程与优化Dynamic flowsheeting and optimization以及用实验数据估计模型参数Estimating model parameters from experimental data。参数估计就是一次特殊的优化目标函数是加权残差平方和、约束是参数取值范围因此初值敏感与不可辨识会复现为初值依赖与局部最优。Q2优化问题必须写清哪三要素少写一个会怎样A决策变量、目标函数与约束。决策变量是允许被搜索的量如进料温度、冷却水流量、温控曲线节点值或参数值目标函数必须是单一标量如能耗最小或残差最小约束含等式约束即模型方程与不等式约束即安全质量与设备上限另需写实施性约束如温度变化率上限。少写决策变量会去优化不该动的量目标函数多目标混写会让求解器无法比较优劣丢掉不等式约束会得到数学最优但物理危险的解例如超温。Q3动态优化与稳态优化在建模与结果形态上有什么本质差别A稳态优化的决策变量是一组常数约束为代数不等式结果是单个最优点用于稳态流程能耗收率优化与设计裕度。动态优化的决策变量是一条轨迹或其若干节点值约束包含必须全程满足的路径约束与终点约束结果是最优曲线与对应状态轨迹用于最优温控曲线与开车停车策略。动态优化有三条硬约束路径约束必须全程满足只在节点检查会漏掉节点之间的超温节点数即曲线分辨率同时决定问题规模与可实施性最优曲线必须能落到现场执行。Q4求解器选项与收敛策略该怎么设为什么只给边界不给经验值A求解器涉及算法选择、收敛判据、可行性容差、迭代上限等大量选项属闭源实现细节名称与默认值在不同版本间存在差异必须以官方文档为准任何数值层面的臆造都会导致读者照抄后得到错误结论。因此只给与实现无关的通行纪律先可行再最优先找一个满足全部约束的解作为起点对量级差异大的变量与约束做无量纲化先简化后加密先用少量决策变量跑通再增加曲线节点宁可多跑几组初值。Q5优化结果为什么要做可信度检查局部最优、约束活性、初值依赖分别怎么查A因为优化器只在你给的模型与约束内寻优不检查就无法区分真的最优与停在局部解。局部最优的查法是用至少三组不同初值复跑并比对目标值差异显著即存在局部最优。约束活性的查法是计算每条不等式约束的间隙间隙为零说明该约束正顶住最优解例如转化率下限恰好活跃意味着结果卡着下限、任何扰动都会使其不合格。初值依赖的查法是记录初值到目标值的对照表并把初值写入结果记录。此外还必须把最优解回代原模型逐条复核不等式约束并做物理合理性复核。本文涉及的语言细节以所购版本随附的官方 gPROMS 文档为准。