混合模型预测控制技术在植物生长室中二氧化碳调节中的应用

《Computers and Electronics in Agriculture》：Hybrid Model Predictive Control for the regulation of carbon dioxide in plant growth chambers

　　吉奥纳塔·奇米尼（Gionata Cimini）|安东尼奥·潘尼科（Antonio Pannico）|马可·加蒂（Marco Gatti）|达尼埃莱·贝尔纳迪尼（Daniele Bernardini）|斯特凡妮亚·德帕斯卡尔（Stefania De Pascale）
ODYS S.r.l.，意大利米兰

**摘要**
精确调节植物生长室内的二氧化碳（CO2）浓度对于确保研究结果的可重复性和生理相关性至关重要。CO2的吸收量会因植物基因型、生长条件、作物密度和物候阶段而显著变化。然而，那些严重依赖作物机制模型的估算和控制方法与植物特征单元（PCU）的目标相悖，因为某些特定作物的模型可能并不可用。此外，在生物再生生命支持系统（BLSS）中，这些方法可能难以应对多种作物、间作、错开收获以及未知的生长阶段。我们提出了一种实时、不依赖于作物特性的方法，通过CO2浓度数据来估算光合作用和呼吸速率。这种方法提高了系统对BLSS典型时变条件的适应性，并支持在没有经过验证的生理模型情况下的操作。所得的速率估算结果可用于诊断算法、监控逻辑和CO2浓度控制器，并为我们的第二个贡献——一种混合模型预测控制（MPC）策略提供了建模基础。该控制器采用混合整数规划方法来处理注射阀的不连续操作范围，并考虑了CO2测量延迟的影响，从而确保在实际运行条件下的质量平衡。我们通过在ESA-MELiSSA框架下实现的PCU中的体内实验验证了这种方法的有效性。

**1. 引言**
生物再生生命支持系统（BLSS）旨在通过减少补给需求和提高可持续性来维持孤立环境（如偏远地区或长期太空任务）中的人类生命。机械系统与生物系统之间的协同作用（受地球生态过程的启发，Wheeler, 2010; Liu et al., 2021）能够生产氧气和食物、去除二氧化碳、净化水以及处理废物（Nelson et al., 2003; Cimini et al., 2025）。BLSS是高度可控的人工生态系统，具有模块化结构，由适当选择的生物组成，能够在封闭的再生循环中执行不同功能（Gitelson and Lisovsky, 2002）。高等植物是BLSS中最关键的组成部分之一，受到了大量研究关注（Zabel et al., 2016）。这是因为高等植物对BLSS的多个需求有重要贡献：它们通过光合作用恢复大气（Sachdeva et al., 2021）、为机组人员提供食物（Carillo et al., 2020）、通过蒸腾作用净化水（Drysdale et al., 2003）、带来心理益处（Szocik et al., 2018），并且可以在一定程度上协助废物处理（Lunn et al., 2017）。这种关键作用与植物与其环境之间的复杂相互作用密切相关。植物在BLSS中的生长和效率取决于许多因素，例如光谱和光照强度、大气成分、空气温度和相对湿度、气流速度以及营养液成分等。因此，应专门开展努力和实验试验，以确定BLSS内各种作物的最佳环境条件（Massa et al., 2016; Nelson et al., 2008）。

植物特征单元（PCU）被设计为研究设施，允许科学家在受控生长条件下研究作物表现，并进行精确的质量平衡测量（Fu et al., 2016）。理想的PCU应具有高度适应性，能够适应不同作物种类，在广泛的气候条件下可靠控制气体交换、植物生理反应和养分吸收动态（Pannico et al., 2022; Peiro et al., 2020; Junker et al., 2015）。其目标是推动硬件设计规范的发展，并确定最佳生长条件和控制方案。在需要控制的众多参数中，二氧化碳（CO2）浓度尤为重要，因为它直接影响光合作用速率和整个系统的生产力（Ainsworth and Long, 2005; Bai et al., 2024）。然而，精确调节CO2浓度面临多个挑战：（i）不同作物的动态响应差异很大；（ii）与其他环境因素（如温度、相对湿度和光照条件）的复杂相互依赖关系会影响植物的光合作用速率；（iii）代谢活动会随着作物的物候阶段而变化。依赖作物动态模型的控制技术可以解决这些挑战。虽然有许多光合作用模型可用（García-Rodríguez et al., 2022），但只有少数模型简单到足以用于控制（Chalabi et al., 2002; Li et al., 2018）。即使那些被认为适用于控制的模型，也主要是为了计算最佳CO2设定点，而且很少有模型被实时应用（Kl?ring et al., 2007）。相反，我们关注的是对确定的、可能随时间变化的CO2轨迹进行稳健和准确的跟踪。在PCU中，精确的CO2控制可以提高实验数据的质量，这些数据可能用于推导复杂的生长模型。在BLSS中，它有助于优化资源利用和与整体优化方案的整合（Ciurans et al., 2022a; Cimini et al., 2025）。Ciurans et al.（2022b）的工作将细胞代谢纳入氧气控制的预测策略中，并提供了仿真结果。

**1.1. 本工作的贡献**
上述所有方法都有一个共同的限制：模型结构和参数都是针对特定目标作物定制的。这意味着它们仅适用于经过广泛测试的作物，而PCU则是专门为研究数据稀缺的作物和/或环境条件而设计的。即使假设所有相关作物在各种环境条件下都有光合作用模型，现实的任务场景可能涉及多种作物的同时种植（间作）、变化的、可能是未知的植物发育阶段以及灵活的收获计划（Tran et al., 2017）。这些因素都会影响现有解决方案的有效性。

我们的第一个贡献是一种实时估算密封室内总CO2交换量的算法，无需对作物类型、数量或生长阶段做出任何假设。该方法使用卡尔曼滤波器进行参数估算，以捕捉主要过程（光合作用和呼吸作用）以及次要效应（如微生物活动、分解或与外部环境的交换）。这一能力本身就是一个重大进步：实时估算使得基于模型的监控、监控逻辑和自适应控制策略成为可能，而无需生理模型。不仅高级控制器可以从这种方法中受益，简单的解决方案（如带有前馈功能的标准比例-积分-微分（PID）调节器）也可以从中受益。通过将速率估算与栽培细节分离，该方法为BLSS中多种作物和生长条件下的更精确气候调节提供了灵活的基础。此外，当次要过程可以忽略或可预测时，该算法可以直接提供光合作用和呼吸速率的可靠计算。

基于我们第一个贡献中的实时速率估算，我们提出了一种模型预测控制（MPC）策略，通过调节注射阀来控制CO2浓度。MPC是这一应用的自然选择，因为它可以明确处理执行限制，并系统地考虑估算的光合作用和呼吸速率作为干扰因素。此外，通过在估算器和预测模型中加入延迟补偿器，我们解决了气体浓度测量延迟这一常见问题（Fitz-Rodríguez et al., 2010; Weiss et al., 2009）。尽管MPC已广泛应用于温室气候控制（Svensen et al., 2024; Morcego et al., 2023a; Yu et al., 2025; Serale et al., 2018）、灌溉（Ghiat et al., 2025）和收获计划（García-Ma?as et al., 2024），但传统系统将CO2注入视为一个昂贵的输入，因为满足最低浓度阈值就足够了。相比之下，PCU需要精确的CO2调节，以实现准确的植物表型和特征实验。这种控制目标的关键差异使我们的工作与众不同，并推动了此处提出的独特MPC公式。

作为MPC公式中的进一步贡献，我们解决了一个据作者所知目前尚未有公开解决方案的限制：流量计通常只在特定操作范围内提供准确测量，而这往往不包括低流量范围。不可靠的流量测量不仅会损害控制性能，还会影响质量平衡的准确性，违背PCU本身的目的。我们的目标是使CO2控制系统始终在这个最佳范围内运行，同时允许控制阀完全关闭。这一要求本质上导致了一个不连续的约束，我们通过混合整数规划方法解决了这个问题（Bemporad and Morari, 1999）。

这两个主要贡献——实时气体交换估算和最佳CO2调节——都通过在那不勒斯费德里科二世大学（Portici, Naples）作物研究实验室的植物特征单元（PCU）中进行的体内实验得到了验证。

**2. 手稿结构**
第2节描述了PCU系统的硬件和软件，并概述了CO2控制问题。第3节开发了光合作用和呼吸活动的估算器，详细说明了如何处理测量延迟及其校准程序。第4节通过体内试验验证了观测器的有效性。第5节介绍了CO2调节的混合MPC公式，并提供了实际实现和校准建议。第6节收集了几项体内实验，以证明控制器的有效性。最后，第7节讨论了未来的研究方向，第8节总结了关键贡献。

**下载：**
下载高分辨率图像（739KB）
下载全尺寸图像

**图1.** 那不勒斯费德里科二世大学作物研究实验室的植物特征单元（PCU）全景视图，该单元是在ESA资助的MELiSSA - PaCMAN项目中实现的。

**2. 系统描述和问题陈述**
在本文中，我们探讨了生长室内CO2控制的关键挑战。尽管这一主题也与农业应用相关，但我们的主要关注点是空间研究。实际上，实验是在那不勒斯费德里科二世大学作物研究实验室（Portici, Naples）的PCU中进行的。该PCU是在ESA - MELiSSA项目“封闭生命支持系统的植物特征单元 - 工程、制造和测试”（PaCMAN）中实现的。图1展示了该系统的全景视图，它由一个2.16立方米的环境控制密封室组成，用于进行准确和广泛的科学研究，为高等植物组件的BLSS集成做准备。有兴趣的读者可以在Pannico et al.（2022）中详细了解该设施。

控制软件部署在一个配备Intel? Core? i3-7100U CPU @2.40 GHz处理器的嵌入式计算平台上。它监控大气和水培子系统，包括超过15个不同的控制回路。这些回路包括空气成分、温度和湿度调节、光照强度、气流、营养液成分和温度以及沟槽液位控制。系统处理超过120个传感器测量值，并控制80多个执行器。控制软件采用多速率设计，最慢的循环周期为2.5秒。在本文中，我们仅关注CO2控制算法，其采样时间为5秒。

CO2浓度κ [ppm]由CAI 700系列NDIR/O2分析仪（California Analytical Instruments, Orange, CA, USA）测量，该分析仪能够测量0至15,000 ppm的CO2。控制器的目标是将κ调节到设定点κ?∈[300,10000] ppm。系统配备了两条平行的CO2注入管线，以适应广泛的CO2输入流量。每条管线都配备了一个质量流量控制器（类型8741，Christian Bürkert GmbH & Co. KG, Ingelfingen, Germany），具有不同的满量程。出于安全考虑，在CO2注入管线分离处安装了一个ON/OFF阀（类型0290，Christian Bürkert GmbH & Co. KG, Ingelfingen, Germany）。系统如图2所示。两条管线设计为独立运行，第一条（管线a）的最大流量为6 NL/h，第二条（管线b）为20 NL/h。出于控制目的，它们被视为一个单一的可操控变量——注入流量，然后根据滞后逻辑分配到相应的管线。设 ηa,k≥0 和 ηb,k≥0 为通往线路 a 和 b 的指令流量，ηk≥0 为控制器发出的流量请求，单位为 NL/h，k 代表离散时间点。我们定义滞后范围 H=[η?b,η?a]。那么线路分配逻辑如下：(1a) 如果 ηk≤η?b，则 ηa,k=ηk；否则 (1b) 如果 ηk≥η?a，则 ηa,k=ηk；否则 ηa,k=0；如果 ηk∈H 且 ηb,k?1=0，则 ηb,k=0；否则 ηb,k=0。η?b 和 η?a 的实际值将在第 5 节中讨论，与控制器设计相关。这种输入映射不影响接下来将讨论的任何控制方法；因此，为了通用性，从现在开始我们主要将 η 称为我们的单一被操纵变量。

下载：下载高分辨率图片（725KB）
下载：下载全尺寸图片
图 2. CO2 注入线的细节，右侧显示了开关阀，左侧显示了两个质量流量控制器。

3. 气体交换率的估计
我们最初感兴趣的是推导 CO2 浓度的一阶摩尔平衡模型。我们假设大气系统的总体积为 vt=4.7m3，包括腔室和循环回路。我们定义 Pa(t) 和 Ta(t) 为循环空气的平均压力 [Pa] 和平均温度 [K]。系统中的气体摩尔数 n(t) 可以通过以下公式计算：(2) n(t)=vtPa(t)RTa(t)，其中 R=8.314J/(mol K)，我们可以推导出测量值 c(t)[mol]=10?6n(t)κ(t) 和设定值 c?(t)[mol]=10?6n(t)κ?(t)。假设气体行为为理想气体，流量注入指令可以转换为 u(t) [mol s-1]：(3) u(t)=Pa,03600Ta,0Rη(t)=k0η(t)，其中 Pa,0=101325Pa 和 Ta,0=273.15K。我们可以定义 ?(t) [mol s-1] 为自由 CO2 交换量，即腔室内所有 CO2 流量的总和，不包括控制注入量。那么摩尔平衡模型为：(4) dc(t)dt=u(t)+?(t)，y(t)=c(t)。显然，使用机械模型估计 ?(t) 非常困难，因为它不仅代表光合作用和呼吸作用，还代表微生物活动、有机物的分解或与外部环境的交换。此外，我们对依赖于作物的模型不感兴趣，因为这可能会削弱 PCU 的目的，特别是在作物变化的情况下。因此，我们采用卡尔曼滤波器在扩展状态空间模型上，其中未知参数被视为具有恒定（随机游走）演变的额外状态，这种方法非常适合于与植物动态相比变化缓慢的参数（Carrassi 和 Vannitsem，2011）。尽管这个假设对于主要由植物生化活动驱动的术语是合理的，但一个关键挑战是光合作用和呼吸作用表现出不同的和切换动态。此外，它们还可能受到环境条件的影响，特别是温度和光谱。我们使用两个不同的缓慢变化参数来估计 ?(t)，这两个参数都表示为 mol s-1：?p(t)∈R，代表光合作用活跃期间的自由 CO2 交换量（称为“模式 p”），以及 ?r(t)∈R，代表呼吸作用活跃期间的交换量（称为“模式 r”）。某些流量可能在两种模式下保持不变，例如由于泄漏到环境中的流量，从而对 ?p(t) 和 ?r(t) 的贡献相同。设 ι(t)∈R [μmol m-2 s-1] 代表腔室内的光照强度，ι? 代表光补偿点，(4) 可以转换为模式依赖的连续时间模型：(5) dc(t)dt=u(t)+?p(t) 如果 ι(t)>ι?；否则 u(t)+?r(t) 如果 ι(t)≤ι?。

为了完全描述系统的混合行为，我们考虑 (6)??ψp(t)?d?p(t)dt,ψr(t)?d?r(t)dt 作为需要估计的两个具有恒定动态的参数，即自由交换量的导数，两者都以 mol s-2 表示。具体来说，选定的模型是一个切换库存系统，其中一个输入被控制（u(t)），而其他输入（?p(t) 或 ?r(t)）作为积分布朗运动演变。这种选择的动机是 PCU 中通常观察到的注入和植物动态之间的明显分离。例如，在本文分析的设施中，从环境条件（约 400 ppm）到典型的 1000 ppm 目标的完全充电需要不到 10 分钟，此时线路 b 固定为 18 NL/h。因此，用于设定点跟踪的增量调整发生在几秒钟的范围内。相比之下，植物吸收的微小变化在更长的时间范围内发生——通常是几小时或更长时间——因此可以认为它们是时间解耦的。我们引入以下扩展状态向量：(7)x(t)=c(t)?p(t)ψp(t)?r(t)ψr(t)? 和简化状态空间：(8a)x(p)(t)=c(t)?p(t)ψp(t)?(8b)x(r)(t)=c(t)?r(t)ψr(t)?，它们分别是 (5) 中两种模式相关状态的分割。假设采样时间为 ts，我们可以定义模式 p 的离散时间状态空间：(9a)xk+1(p)=A(p)xk+B(p)uk+ωk(p)(9b)yk=C(p)xk(p)+vk(p) 和模式 r 的状态空间：(10a)xk+1(r)=A(r)xk+B(r)uk+ωk(r)(10b)yk=C(r)xk(r)+vk(r)，其中 A(p), A(r)∈Rne×ne, B(p), B(r)∈Rnu×ne 和 C(p), C(r)∈Rny×ne。从 (8) 我们得到 ne=3, nu=1, ny=1。公式 (9) 中的 ω(p) 和 v(p) 分别代表过程噪声和测量噪声，它们模拟了状态空间方程的不确定性和测量误差。我们假设这些是白噪声，因此 ωk(p)～N(0,Q(p)) 和 vk(p)～N(0,R(p))，其中 R(p)>0, Q(p)≥0。类似的考虑也适用于 (10)，因此 ωk(r)～N(0,Q(r)) 和 vk(r)～N(0,R(r))，其中 R(r)>0, Q(r)≥0。从模式依赖模型 (5) 和状态空间选择 (8) 可以容易地看出 A(p)≡A(r)≡A, B(p)≡B(r)≡B 和 C(p)≡C(r)≡C 成立，其中 (11)Ac=010001000,B=ts00,C=100，且 A=etsAc。

我们可以为 xk(m) 导出两个最优的稳态估计器，其中 m∈{p,r} 是模式标签。我们引入符号 x?k|k(m)，它是给定时间 k 时 x(m) 的最优估计，基于直到时间 k 的测量值。稳态卡尔曼方程为：(12a)x?k|k(m)=x?k|k?1(m)+M(m)(yk?Cx?k|k?1(m))(12b)x?k+1|k(m)=Ax?k|k(m)+Buk，其中稳态增益 M(m)∈Rne×ny 可以离线计算为：(13)M(m)=P(m)C?CP(m)C?+R(m)?1，P(m) 是代数 Riccati 方程的唯一正定解：(14)P(m)=AP(m)A?+Q(m)?AP(m)C?CP(m)C?+R(m)?1CP(m)A?。利用 (12)，我们可以仅为给定模式的相关状态计算估计值，这些状态由 ιk 的值决定（见 (5)）。在正常操作中，从模式 p 切换到模式 r 对应于白天/黑夜的转换，反之亦然。每种模式保持活跃的时间足够长，以便在新相关状态的最后一个可用估计值在切换发生之前变得过时。例如，设 k1 和 k2 分别是光照到黑暗阶段和黑暗到光照阶段的时刻，且 k2>k1。由于植物在黑暗阶段生长，当灯光重新打开时消耗的 CO2 量将超过 ?p 的最后一个可用估计值，使得 ??p,k1 不是 ??p,k2 的准确初始猜测。唯一在模式变化过程中保持一致的状态是 ck，它被两个简化状态空间模型 (8) 共享。因此，我们的想法是在切换时刻通过开环演化的自由 CO2 交换提供可靠的初始猜测，这没有被活动模式捕获。这就是在扩展状态空间中引入流量导数的原因。我们在第 3.2 节提供了估计器的完整公式，在此之前我们详细介绍了处理测量延迟的方法。

3.1. 气体测量延迟
众所周知，腔室响应时间和气体-壁分离显著影响气体测量，通常会引入不可忽视的延迟（Deming 等人，2019；Weiss 等人，2009）。为了考虑这种效应，我们提出了一个考虑延迟的 (22) 公式，以允许更精确和可靠的估计。设 τc 为 CO2 测量的延迟，使得 y(t)=c(t?τc)。我们的目标是确定延迟步数 τ=?τc/ts?∈Z，其中 ??? 是向最接近整数的四舍五入运算符，Z 是正整数集。我们专注于离线计算 τ 的程序。互相关最大化是最广泛采用的方法之一。考虑 (4)，由于 ?(t) 的带宽远低于 u(t)，它对时间对齐的贡献很小。因此，即使无法进行空腔室的专门测试，从生活实验中记录的数据仍然适用。为了计算平滑导数 dy(t)/dt 的互相关，我们使用基于局部多项式拟合的微分滤波器（即 Savitzky–Golay 滤波器）。假设 y(t) 以等间隔采样，采样时间为 ts。对于微分滤波器，考虑一个长度为 w?γ=2wγ+1 和多项式度数为 pγ 的窗口，其中 pγ≤w?γ?1。设 V∈Rwγ×(pγ+1) 表示设计矩阵，计算如下：(15)Vi,j=(i?wγ?1)j?1,i=1,…,w?γ,j=1,…,pγ+1，并定义 G∈Rw?γ×(pγ+1) 为微分系数矩阵，由 G=V(V?V)?1 给出。我们可以通过以下卷积计算 yk 的平滑一阶导数，表示为 y?k′∈R：(16)y?k′=1ts∑i=?wγwγGi+wγ+1,2yk+i。由于 G 的第 i 列提供了 i?1 导数系数。对于一般的时间序列 X∈RN，其中 X=[x1,…,xN]?，我们引入归一化运算符 ζ：R→R，它将 xi, i=1,…,N 映射为其对应的 z 分数标准化值 x?i，即：(17)x?i=ζ(xi,X)=xi?μ(X)σ(X)，其中 μ,σ:RN→R 分别是时间的均值和标准差运算符，使得：(18)μ(X)=1N∑i=1Nxi,σ(X)=1N∑i=1N(xi?μ(X)2。使用之前介绍的仪器，我们可以确定延迟步数 τ。考虑数据集 {(Yj,Uj)，j=1,…,nt}，其中 Uj, Yj∈Rmj 分别代表 ui 和 yi 的时间序列，i=1,…,mj。对 (Yj,Uj) 的配对可能代表同一实验的部分或来自不同的实验。自然地，时间序列 Y?j′ 是通过将方程 (16) 应用于 Yj 而得到的。给定 D∈Z 为候选延迟值集，τ 的最优值通过解决以下互相关最大化问题确定：(19)τ=argmaxd∈D1nt∑j=1nt1mj?d∑i=d+1mjζ(y?i′,Y?′)ζ(ui?d,U)。这种延迟不一定随时间保持不变。例如，环境条件的重大变化，特别是空气流量，可能会产生时变行为。因此，可以假设 τk=τ+τ?k，其中 τ∈Z 是使用 (19) 离线计算的时间不变分量，τ?k∈Z 是在线学习的时变校正。一种可能的方法是使用方程 (19) 在线计算更新的延迟估计 τg,k。然后，为了减轻抖动，可以通过对 τk?τg,k 的差值进行指数平滑来获得 τ?k。

3.2. 实际考虑和实施改进
为了改善数值条件并使积分器链的三个组成部分在幅度上可比，我们引入了一个特征时间尺度 β。具体来说，我们定义缩放矩阵：(20)S=1000β000β2，β?1，并通过 x(m)=Ss(m) 定义缩放状态 s(m)，使得 s(m)=cβ?mβ2ψm?。这里，β 可以解释为“时间单位”，即 c 有意义变化的特征时间范围。等效地，用无量纲时间 ν=t/β 表示，缩放状态 β?(m) 和 β2ψ(m) 分别对应于 c 关于 ν 的一阶和二阶导数分量。如果我们想用缩放坐标表示离散时间模型，很容易验证 CS≡C 和 S?1B≡B，因此 yk=Csk(m)+vk(m) 显然成立，而离散时间状态演变为 sk+1(m)=ASsk(m)+Buk+ωS(m)，其中 AS=S?1AS，ωS(m)～N(0,QS(m)) 和 QS(m)=S?1QS??。因此，我们引入 PS(m)，它是缩放坐标下代数 Riccati 方程的唯一正定解：(21)PS(m)=ASPS(m)AS?+QS(m)?ASPS(m)C?CPS(m)C?+R(m)?1CPS(m)AS?。通过 (12)，我们可以仅为给定模式的相关状态计算估计值，这些状态由 ιk 的值决定（见 (5)）。在正常操作中，从模式 p 切换到模式 r 对应于白天/黑夜的转换，反之亦然。每种模式保持活跃的时间足够长，以便新相关状态的最后一个可用估计在切换发生之前变得过时。例如，设 k1 和 k2 分别是光照到黑暗阶段和黑暗到光照阶段的时刻，且 k2>k1。由于植物在黑暗阶段生长，当灯光重新打开时消耗的 CO2 量将超过 ?p 的最后一个可用估计值，使得 ??p,k1 不是 ??p,k2 的准确初始猜测。唯一在模式变化过程中保持一致的状态是 ck，它被两个简化状态空间模型 (8) 共享。因此，我们的想法是在切换时刻通过开环演化的自由 CO2 交换提供可靠的初始猜测，这没有被活动模式捕获。这就是在扩展状态空间中引入流量导数的原因。我们在第 3.2 节提供了估计器的完整公式，在此之前我们详细介绍了处理测量延迟的方法。

3.1. 气体测量延迟
众所周知，腔室响应时间和气体-壁分离显著影响气体测量，通常会引入不可忽视的延迟（Deming 等人，2019；Weiss 等人，2009）。为了考虑这种效应，我们提出了一个考虑延迟的 (22) 公式，以允许更精确和可靠的估计。设 τc 为 CO2 测量的延迟，使得 y(t)=c(t?τc)。我们的目标是确定延迟步数 τ=?τc/ts?∈Z，其中 ??? 是向最接近整数的四舍五入运算符，Z 是正整数集。我们专注于离线计算 τ 的程序。互相关最大化是最广泛采用的方法之一。考虑 (4)，由于 ?(t) 的带宽远低于 u(t)，它对时间对齐的贡献很小。因此，即使无法进行空腔室的专门测试，从生活实验中记录的数据仍然适用。为了计算平滑导数 dy(t)/dt 的互相关，我们使用基于局部多项式拟合的微分滤波器（即 Savitzky–Golay 滤波器）。假设 y(t) 以等间隔采样，采样时间为 ts。对于微分滤波器，考虑一个长度为 w?γ=2wγ+1 和多项式度数为 pγ 的窗口，其中 pγ≤w?γ?1。设 V∈Rwγ×(pγ+1) 表示设计矩阵，计算如下：(15)Vi,j=(i?wγ?1)j?1,i=1,…,w?γ,j=1,…,pγ+1，并定义 G∈Rw?γ×(pγ+1) 为微分系数矩阵，由 G=V(V?V)?1 给出。我们可以通过以下卷积计算 yk 的平滑一阶导数，表示为 y?k′∈R：(16)y?k′=1ts∑i=?wγwγGi+wγ+1,2yk+i。由于 G 的第 i 列提供了 i?1 导数系数。对于一般的 time-series X∈RN，其中 X=[x1,…,xN]?，我们引入归一化运算符 ζ：R→R，它将 xi, i=1,…,N 映射为其对应的 z 分数标准化值 x?i，即：(17)x?i=ζ(xi,X)=xi?μ(X)σ(X)，其中 μ,σ:RN→R 分别是时间的均值和标准差运算符，使得：(18)μ(X)=1N∑i=1Nxi,σ(X)=1N∑i=1N(xi?μ(X)2。使用之前介绍的仪器，我们可以确定延迟步数 τ。考虑数据集 {(Yj,Uj)，j=1,…,nt}，其中 Uj, Yj∈Rmj 分别代表 ui 和 yi 的时间序列，i=1,…,mj。对 (Yj,Uj) 的配对可能代表同一实验的部分或来自不同的实验。自然地，时间序列 Y?j′ 是通过将方程 (16) 应用于 Yj 而得到的。给定 D∈Z 为候选延迟值集，τ 的最优值通过解决以下互相关最大化问题确定：(19)τ=argmaxd∈D1nt∑j=1nt1mj?d∑i=d+1mjζ(y?i′,Y?′)ζ(ui?d,U)。这种延迟不一定随时间保持不变。例如，环境条件的重大变化，特别是空气流量，可能会产生时变行为。因此，可以假设 τk=τ+τ?k，其中 τ∈Z 是使用 (19) 离线计算的时间不变分量，τ?k∈Z 是在线学习的时变校正。一种可能的方法是使用方程 (19) 在线计算更新的延迟估计 τg,k。然后，为了减轻抖动，可以通过对 τk?τg,k 的差值进行指数平滑来获得 τ?k。

3.2. 实际考虑和实施改进
为了改善数值条件并使积分器链的三个组成部分在幅度上可比，我们引入了一个特征时间尺度 β。具体来说，我们定义缩放矩阵：(20)S=1000β000β2，β?1，并通过 x(m)=Ss(m) 定义缩放状态 s(m)，使得 s(m)=cβ?mβ2ψm?。这里，β 可以解释为“时间单位”，即 c 有意义变化的特征时间范围。等效地，用无量纲时间 ν=t/β 表示，缩放状态 β?(m) 和 β2ψ(m) 分别对应于 c 关于 ν 的一阶和二阶导数分量。如果我们想用缩放坐标表示离散时间模型，很容易验证 CS≡C 和 S?1B≡B，因此 yk=Csk(m)+vk(m) 显然成立，而离散时间状态演变为 sk+1(m)=ASsk(m)+Buk+ωS(m)，其中 AS=S?1AS，ωS(m)～N(0,QS(m)) 和 QS(m)=S?1QS??。因此，我们引入 PS(m)，它是缩放坐标下代数 Riccati 方程的唯一正定解：(21)PS(m)=ASPS(m)AS?+QS(m)?ASPS(m)C?CPS(m)C?+R(m)?1CPS(m)AS?。通过将 P(m) 替换为 PS(m) 在 (13) 中获得缩放的稳态增益 MS(m)。设延迟时间步长 ?=k?τk。符号 x??|k（或在缩放坐标的情况下为 s??|k）表示在时间步长 ? 给定所有测量信息时 x（或 s）的估计值。给定 (12) 和非相关状态的开环假设，我们得到以下延迟估计器：(22)ifι?>ι?，s??|k(p)=s??|k?1(p)+MS(p)(yk?Cs??|k?1(p))s??+1|k(p)=ASs??|k(p)+Bu???r,?+1|k=??r,?|k?1+tsψ?r,?|k?1ψ?r,?+1|k=ψr,?|k?1；ifι?≤ι?，s??|k(r)=s??|k?1(r)+MS(r)(yk?Cs??|k?1(r))s??+1|k(r)=ASs??|k(r)+Bu???p,?+1|k=??p,?|k?1+tsψ?p,?|k?1ψ?p,?|k?1。注意，当 ι?>ι? 时，??r,?|k≡??r,?|k?1 和 ψ?r,?|k≡ψ?r,?|k?1，因为在这种情况下状态演算是开环的。类似地，当 ι?≤ι? 时，??p,?|k≡??p,?|k?1 和 ψ?p,?|k≡ψ?p,?|k?1。显然，如果没有测量延迟，则 ?≡k（τk≡0）。在 (22) 中的切换卡尔曼滤波器提供了完整状态向量 (23)s??|k=c??|kβ??p,?|kβ2ψ?p,?|kβ2ψ?r,?|k? 的估计，从中可以容易地导出 (24)x??|k=c??|k??p,?|kψ?p,?|k??r,?|kψ?r,?|k?。请注意，标准生长室协议通常使用方波辐照度曲线，其中光照强度瞬间切换（Poorter 等人，2016）。在这些条件下，这也适用于目标设施中使用的研究协议，ι? 的设置很直接。如果使用更复杂的曲线，例如正弦曲线（Stamford 等人，2024），则需要相应地修改切换函数，可能需要考虑瞬态行为。然而，这超出了本文的范围，因为底层的估计器和控制器保持不变。我们现在提供了对 (22) 的潜在简化方法的见解，这些简化可能对在极低成本单元上的实现有益。光照和黑暗阶段之间的环境变化可能对 CO2 传感器和过程动态的不确定性影响可以忽略不计。在这种假设下，我们可以施加 (25a)?QS(p)≡QS(r)?QS(25b)?R(p)≡R(r)?R，由此，根据 (13) 和 (14)，可以很容易地得出 (14)，因此公式简化为存储单个增益矩阵。基本控制设置通常针对图 6（顶部图表）中描绘的锯齿形 CO2 曲线。在植物的光照阶段，CO2 水平保持在规定的设定点。在黑暗阶段，植物通过呼吸作用产生 CO2，控制系统保持不活跃，CO2 水平上升。随着光照阶段的恢复，光合作用重新开始，CO2 减少，直到再次达到设定点，从而触发控制器重新激活。在这种简化的操作中，黑暗阶段不计划进行控制干预。因此，如果仅为了支持控制动作而估计 ?(t)，则估计可以限制在 ?p(t) 上。此外，由于在黑暗阶段之后二氧化碳水平升高，估计器可以从整个消耗阶段中受益，以收敛到实际值，从而减少了对开环近似的需要。虽然这些简化不可避免地降低了估计器的可靠性、灵活性和鲁棒性，但当控制要求简单时，它们代表了一个有效的替代方案，具有实现更轻量级的优势。

3.3. 估计器校准
本节概述了一种系统化的、特定于应用的卡尔曼滤波器参数调整程序。在整个过程中，我们假设（25）成立，然而该方法可以直接扩展到使用不同的协方差矩阵来表示光照和黑暗时期的情况。对于植物生长实验中典型的锯齿形 c 曲线（在第 3.2 节中解释），我们引入了第 j 天消耗期（光照开始时）和积累期（黑暗时期）的二氧化碳变化率的平均值，分别表示为 ??p,j 和 ??r,j，定义为：
(26a) ??p,j = 1/tp,j + ∫tp,j?tp,j? c?(ν)dν
(26b) ??r,j = 1/tr,j + ∫tr,j?tr,j? c?(ν)dν
其中 [tp,j?,tp,j+] 表示消耗时间间隔，[tr,j?,tr,j+] 表示积累时间间隔。显然，这些间隔是连续发生的，并且 tr,j+≡tp,j+1?。在实践中，??p,j 和 ??r,j 代表构成开环曲线 c 的斜率的平均值，分别以中点 ?tp,j?(tp,j++tp,j?)/2 和 ?tr,j?(tr,j++tr,j?)/2 为中心。为了获得稳健的导数，我们对整个间隔应用 Savitzky–Golay 滤波器，类似于（19）。由于这样的估计本质上是反因果的，它们不适用于实时使用。尽管如此，它们非常适合于观察器调整和离线验证。

可以从传感器数据表中获得更原理性的测量噪声协方差 R 的估计。然而，制造商通常指定最坏情况界限或峰对峰值噪声水平。这些值在设计上是保守的，很少能代表最佳调整所需的随机特性。为了获得基于经验的测量噪声特性描述，我们使用了空腔在稳态条件下的数据，从而确保测量的信号 ck 反映的是传感器噪声而不是物理动态。大约使用了 10 小时的连续数据，从中计算了测量波动的样本标准差 σ(c)。假设噪声为零均值且平稳，然后将测量噪声协方差设置为 R=σ(c)2=8?10??，这提供了一个针对传感系统实际运行情况的统计上有意义的估计。

在调整协方差矩阵 QS 时，我们利用了腔内二氧化碳快速动态与植物较慢生物过程之间的明显时间尺度分离。由于 ψ(m) 的时间演化被建模为随机游走，其变异性必须通过过程噪声明确捕获。同时，我们要求注入的不确定性必须一致地传播到更高阶状态，确保沿积分器链的连贯传输。为了保留强制所需时间尺度分离的自由度，我们引入了一个直接作用于 c 的额外噪声分量。我们定义连续时间谱密度 Qc=diag(10q1,10qrq1)，其中 ?q?[q1,qr]∈R2 是校准向量，且 qr<0。相应的离散时间协方差矩阵通过（27）获得：
?QS = ∫?^t seAνGQcG?eA?νdν，G?100001，
其中选择矩阵 G 确保正确的噪声分配。让 Df 表示用于校准的体内实验中收集的 Nf 个样本的数据集，包括测量和控制输入。为了简化符号，我们用 x?(qf) 表示由卡尔曼滤波器产生的通用状态 x 的估计，其稳态增益由 qf 唯一确定；R 如前一段所述是固定的。对于任何时间序列 v，我们定义增量运算符 δ(vk)=vk?vk?1。我们建议根据在 Df 上计算的三个性能指标来调整 q：
g1：与二氧化碳浓度估计相关的创新误差；
g2：与平均反因果生物速率 ??m 的差异，如（26）中计算的那样，覆盖 D 个锯齿形曲线（D 通常对应于实验天数）；
g3：惩罚光合作用和呼吸活动导数快速变化的平滑度指标。

这些性能指标在以下公式中形式化：
(28a) g1(qf) = 1/Nf?τ ∑k=τ^Nf(yk?c??|k?1(qf))2
(28b) g2(qf) = 1/2D ∑j=1^D(??p,j???r,j?|j?(qf))2 + (??r,j???r,j?|j?(qf))2
(28c) g3(q) = 1/2(Nf?τ?1) ∑k=τ^Nf?1δ(??p,?|k(qf))2 + δ(??r,?|k(qf))2。
该方法背后的理由很简单：观察器必须准确预测一步之后的二氧化碳浓度，同时匹配潜在的呼吸和光合作用速率，偏好平滑的随机游走行为，以抑制估计活动中生物学上不合理的波动。我们的目标是解决多目标优化问题：
(29) min qfg(qf), g(qf) = [g1(qf), g2(qf), g3(qf)]?。
由于目标是竞争的，解决方案在帕累托意义上是有意义的。在实践中，一种有效且透明的方法是在一组候选 (q1,qr) 对上评估 g(qf) 并保留非支配点。图 3 显示了在几个帕累托最优选择 (q1,qr) 处评估的向量值函数 g(q)。采用的调整是 q1=?2, qr=?5，这转化为谱密度 Qc=diag(10?2,10??)。值得强调的是，这里描述的校准程序完全可以在回放数据上执行，无需实现估计器，并且独立于所使用的控制器。

下载：下载高分辨率图像（262KB）
下载：下载全尺寸图像

图 3. 对于定义卡尔曼滤波器增益 MS 的帕累托最优选择 (q1,qr)，g(qf) 的标准化值（见公式（28）。
采用的实现选择是 q1=?2, qr=?5。

表 1. 生命测试设置列表。

ID
作物
持续时间 [d]
昼夜 [h]
空气温度（昼夜）[°C]
空气相对湿度（昼夜）[%]
二氧化碳设定点 [ppm]
实验 1
生菜（Lactuca sativa L.）
28
16/8
26/20
50/70
1000
实验 2
羽衣甘蓝（Brassica oleracea L. var. acephala）
28
16/8
26/20
50/70
1000

4. 估计器的实验验证
我们在第 2 节描述的 PCU 进行的两个实验中测试了（22）的估计性能。这些实验属于 PaCMAN 项目的活动范围。设置收集在表 1 中，图 5 显示了实验结束时两种不同作物的状态。

首先，为了离线估计 τ，我们从之前的实验数据中提取了 nt=30 对输入-输出对 (Yj,Uj)。成本函数（19）的不同值 d∈{0,…,14} 收集在图 4 中，最大值 τ=7 被突出显示。在我们的设置中，环境条件不预期会导致与离线估计的延迟有显著偏差；因此，我们假设 τk=τ。

两个实验的估计算法结果如图 6 所示。对于两者，设定点 κ?=1000ppm 在实验过程中保持不变。由于腔内没有配备从大气中去除二氧化碳的仪器，κ? 应被视为在光照阶段需要维持的水平，从而导致顶部图表的典型锯齿形曲线。值得注意的是，在这种情况下，可以考虑仅依赖 ??p 和 ??r。然而，这只能提供有限的效用，因为它产生的估计是偶发的（每天一次）且反因果的，不适合动态环境条件或收获。即使在没有不可预见事件的情况下，植物在光照阶段的二氧化碳消耗也可能显著变化（如图 6 所示，超过 30%），因此这种简单的导数方法不适合辅助控制任务或用于实时诊断算法。尽管如此，它仍然有助于验证实时估计算法（22）的性能。这正是我们在图 6 的第二行图表中所使用的，其中 ??p 和 ??r 代表 ? 的实时估计，与根据（26）和后续描述离线计算的 ??p 和 ??r 进行比较。估计参数的演变是平滑的，与实验过程中的离线反因果值 ??p 和 ??r 相匹配。由于后者只是几小时窗口内的平均（过滤）导数，因此与实时估计的小偏差是预期的。估计误差 |c??c??|k| 平均低于 10??，因此图中没有报告。

下载：下载高分辨率图像（137KB）
下载：下载全尺寸图像

图 4. 在一组 nt=30 个输入-输出对 (Yj,Uj) 上计算的交叉相关值。显示了集合 D=[1,…,14] 的成本函数（19）的值，并突出显示了对应于 τ 选择的最大值。图 6 显示，延迟估计器（22）在两个实验中提供了可靠的 ??p 和 ??r 估计，正如我们在下一节中讨论的，这两个实验使用了两种非常不同的控制器。这些实时的、与作物无关的估计构成了我们的第一个贡献：它们仅依赖于二氧化碳浓度数据，可以支持系统监控、增强控制性能和异常诊断。与图 6 的顶部图表相关的控制方面将在下一节中讨论。

下载：下载高分辨率图像（526KB）
下载：下载全尺寸图像

(a). 实验 1：生菜（Lactuca sativa L.）
下载：下载高分辨率图像（480KB）
下载：下载全尺寸图像

(b). 实验 2：羽衣甘蓝（Brassica oleracea L. var. acephala）。

图 5. 实验结束后立即在腔内观察到的作物情况，以及收获前。

下载：下载高分辨率图像（420KB）
下载：下载全尺寸图像

(a). 实验 1。
下载：下载高分辨率图像（441KB）
下载：下载全尺寸图像

(b). 实验 2。

图 6. 在那不勒斯费德里科二世大学 ESA-MELiSSA PCU 进行的两个实验中状态估计器（22）的性能评估（左图为第一个实验，右图为第二个实验），见表 1。顶部图表显示了 κ 的演变，中间图表显示了光合作用和呼吸活动的实时估计，底部图表完成了卡尔曼滤波器的完整状态向量（24）和恒定动态参数。

5. 混合 MPC 用于精确的二氧化碳控制
准确的实时气体交换率估计几乎可以改进任何二氧化碳调节策略。即使是无模型的调节器也可以通过前馈动作利用这些估计，如结果后面所示。然而，在 PCU 中，当有足够的计算资源时，MPC 是最自然的选择。首先，MPC 可以将估计的植物活动项 ??p,k 和 ??r,k 直接作为测量干扰纳入预测模型中，当植物驱动的动态占主导地位并且在不同作物和生长阶段变化时，这可以提高性能。其次，MPC 明确地强制执行输入、状态和输出约束，使其非常适合执行器限制，这在库存型系统（如气体注入）中经常受到强调。第三，MPC 自然地适应了第 3.1 节讨论的测量延迟，并且最后，它允许将控制目标表述为将二氧化碳维持在规定的范围内，而不是跟踪单一设定点。这在 BLSS 操作中可能是有利的，例如更有效地协调多个隔间并提高整体系统性能。

操作上，MPC 将控制问题重新表述为一个优化问题，在每个时间步骤中，计算一个最优的操作序列（满足输入、状态和输出约束）。然后只将序列中的第一个操作应用于系统。BLSS 中的复杂控制任务，如隔间的协调（Ciurans 等，2022a）或气候调节（Blasco 等，2007，Morcego 等，2023a，Yu 等，2025）最近受益于基于模型的技术的使用。

可靠地测量控制算法注入的二氧化碳流量对于准确的质量平衡计算至关重要。然而，流量计在阀门位置较窄时往往不准确，因此尽可能避免这种情况是可取的。为了解决这个实际约束，我们在 MPC 公式中对操纵变量施加了一个不相交的可行集，确保阀门要么完全关闭（零），要么在保证测量准确性的最小阈值以上运行。考虑操纵变量 u，我们想施加不相交的约束（30）uk∈{0}∪[umin,umax]，其中 umin=k0ηmin，umax=k0ηmax，并且 0< />ι?
??r,?+i|k
if i?+i≤ι?
其中 i=2,…,τk，初始条件 ??p,?+1|k, ψ?p,?+1|k, ??r,?+1|k, ψ?r,?+1|k 和 c??+1|k 从（22）计算得出，导数 ??p,?+1|k, ??r,?+1|k 保持不变（见（31b），（31c））。迭代31次，每次τk步，得到以下状态向量估计：(32)x?k|k=c?k|k??p,k|kψ?p,k|k??r,k|kψ?r,k|k，用作MPC预测的初始条件。在(32)中，ψ?p,k|k≡ψ?p,?+1|k和ψ?r,k|k≡ψ?r,?+1|k显然成立。设Wδ∈Rnu×nu和Wy∈Rny×ny为输入率和输出跟踪的正权重矩阵。在我们的案例中，Wδ和Wy都是标量，因为nu=ny=1。为了简化符号，让索引j表示时间步k的第i个预测k+i，其中j,i=0,…,p?1，p是预测步数。我们定义第j阶段成本为以下二次函数：(33)J(uj,ej+1)=‖Wδδ(uj)‖22+‖Wyej+1‖22，其中?ej?c?j?cj是跟踪误差，p是预测范围。通过(33)，我们希望同时惩罚偏离CO2设定点和阀门开度的频繁变化。所提出的MPC旨在解决以下优化问题：(34a)min{cj+1,uj}j=0p?1∑j=0p?1J(δj,ej+1)(34b)s.t.cj+1=cj+ts(uj+??j|k)(34c)uj∈{0}∪[umin,umax]，其中??j|k是预期对状态动态的干扰。我们通过将优化问题重新表述为混合整数规划（MIP）的一个子类——混合二进制二次规划（MB-QP）问题来解决(34)的非凸性。为此，我们引入p个辅助优化变量ρj∈{0,1}，并将(34c)重新表述为：(35)ρjumin≤uj≤ρjumax,j=0,…,p?1。为了降低(34)的维度，我们应用状态消除技术，也称为压缩（Bemporad和Cimini，2023）。这种技术递归地将预测步骤中的状态序列c表示为控制输入u和测量干扰?的函数。因此，可以消除等式约束，允许直接在定义为z=uρ的简化决策变量集上执行优化，其中u=u0…up?1且ρ=ρ0…ρp?1。设n=2p，并考虑向量θ∈Rn+2，收集问题(34)的所有时变参数，定义为：(36)θk=uk?1ck??k|k…??k+p?1|kc?k+1…c?k+p′。在每个时间步，我们必须解决以下压缩的参数MB-QP问题：(37a)?minz?[u,ρ]0.5z′Hz+(Fθk)′z(37c)s.t.Az≤bρ∈{0,1}，其中H∈Rn×n是半正定的Hessian矩阵，F∈Rn×(n+2)，A∈Rn×n且b∈Rb。不等式约束Az≤b的定义很容易从(36)得出：(38a)A=Ip?umaxIp?IpuminIp(38b)b=0…0，其中Ip∈Rp×p表示单位矩阵。5.1. 控制器校准和实现细节考虑目标函数(34a)，我们首先描述一个适当校准参数Wy、Wδ和p的程序。优化问题的解对成本缩放是不变的，因此，由于nu=ny=1，我们可以将权重调整限制为标量?w?Wδ/Wy∈R。由于(34b)具有无限的直流增益，预测范围不能根据稳定时间来调整，而必须(i)满足预览要求，(ii)作为w之外的调节旋钮，用于闭环带宽。因此，我们提出了一种类似于第3.3节中的方法来适当选择?qc?[w,p]。设Dc表示在MPC控制器(34)计算的作用下，对(9)–(10)的闭环演变收集的Nc个样本的数据集。我们的目标是解决argminqch(qc)，其中解决方案在帕累托意义上是有意义的。向量值函数定义为h(qc)=[h1(qf),h2(qc)]?，其中：(39a)h1(qc)=1Nc∑k=0Nc?1ek(qc)2(39b)h2(qc)=1Nc∑k=0Nc?1δ(uk(qc))2，我们用ek(qc)、δ(uk(qc))表示给定qc选择下的闭环演变的跟踪误差和输入增量。每次评估h(qc)都需要在多周内收集新的闭环数据，以捕捉植物吸收动态的效果。因此，与可以在回放模式下使用真实数据进行的观察器调整不同，控制器调整依赖于模拟来生成所需的闭环轨迹。由于?(m)是测量到的干扰，因此不受闭环的影响，我们建议使用真实数据或其简化模型来表示其轮廓。在这项工作中，我们确定了两个单独的Sigmoid函数，sp(t)和sr(t)，分别拟合实验1中?(p)(t)和?(r)(t)的轮廓。然后从三个闭环模拟中构建数据集Dc，在这些模拟中，每个Sigmoid轮廓的终端值被缩放到其标称值的1、2和3倍。这个程序探索了一组具有代表性和挑战性的干扰演变，同时保持了原始轮廓的定性形状。预测范围的适当调查范围应与注入/混合动态兼容，因此第3.1节中的延迟分析是识别正确候选范围的良好基础。更长的（生物）时间框架范围内的预测范围不会带来控制改进，因为动态在时间上很大程度上是解耦的（见第3节）。图7显示了不同qc对选择的指标(39)。由于更大的预测范围对应于增加的计算复杂性，如果性能没有显著提高，则更倾向于选择较低的值，因此在这项工作中我们选择了p=6和w=10?4。问题(37)在数学上等同于(34)，优化器z可以用任何MIQP求解器计算。这些方法通过分支定界、切割平面或混合启发式（Achterberg等人，2005）来处理连续和离散变量的复杂性。在这项工作中，我们使用“ODYS QP求解器”（Cimini等人，2017），以下简称ODYSQP，其最新版本支持连续和二进制变量。它是一个高性能且可靠的求解器，专门用于实时嵌入式优化中出现的问题。为了与最先进的求解器进行吞吐量基准测试，并且通常用于MPC应用，我们建议感兴趣的读者参考Cimini和Bemporad（2025）。设置p=6，相当于30秒的预测范围，总内存占用为4.94 kB，包括问题数据和ODYSQP内部分配。为了测试在第2节描述的控制单元上实时解决(37)的计算可行性，我们从实验1的数据集中收集了θk集合，大约有250万个问题，包括6个连续优化变量、6个离散优化变量和12个约束。使用ODYSQP解决(37)的最坏情况周转时间为233微秒，这对整个PCU控制软件的影响可以忽略不计。下载：下载高分辨率图像（196KB）下载：下载全尺寸图像图7. 对于定义卡尔曼滤波器增益MS的调优对(q1,qr)的帕累托最优选择的h(qc)的标准化值（见方程(28)）。采用的实现选择是q1=?2，qr=?5.6。控制器的实验验证我们通过比较表1中两个实验的结果来研究MPC(34)的性能，尽管这两个实验运行了相同的估计器(22)，但实现了两种不同的控制器。在第一个实验中，CO2水平由一个相当标准的PI控制器调节，该控制器具有饱和控制输入和通过反饱和方法获得的反饱和动作。设ek表示时间步k的跟踪误差。饱和控制输入uk通过以下操作序列计算：(40a)u?k=Kpek+Ik(40b)uk=max(min(u?k,umax)，umin)(40c)Ik+1=Ik+tsKiek+Ka(uk?u?k)，其中Ik是累积积分动作（I0=0），Kp、Ki和Ka分别是比例、积分和反饱和动作的正标量权重，u?k表示未饱和的控制动作。对于CO2调节的标准实现(40)，如Pannico等人（2022）中详细描述的，我们可以使用以下跟踪误差：(41)ek=c?k?yk。如前所述，如果ek按照(41)计算，根据测量延迟对yk和植物气体交换率的影响，(40)的性能可能会很差。但是，估计器(22)也可以用来增强比MPC更基本的控制器。考虑另一种跟踪误差(42)ek=c?k?c?k|k，其中c?k|k是首先使用(22)计算得到的x??|k。这个估计值作为(31)的初始条件，进而得到(32)中的x?k|k。在本文报告的两个实验中，我们测试了上述所有讨论的CO2调节控制提案，具体来说：•在实验1中，从第1天到第18天的一半，使用PI(40)和跟踪误差(41)；•在实验1中，从第18天的一半到实验结束，使用PI(40)和跟踪误差(42)。我们将这种替代方案称为高级PI；•实验2完全使用第5节中的混合MPC。我们回想一下，在实验中施加了一个常数κ?=1000ppm。CO2的总体演变可以在图6的顶部图表中观察到，我们感兴趣的是比较跟踪误差|κ??κk|。图8收集了这种绝对误差的统计洞察，通过箱形图汇总，范围从第25百分位数到第75百分位数。从用于计算分布的数据中，我们忽略了与黑暗小时和随之而来的光照不足期相关的数据间隔，因为这些偏差与使用的控制器无关，无论在那些间隔期间控制器都是不活跃的。两个额外的事件导致CO2水平的预期偏差。第一个事件涉及光照辐射。每小时，安装在PCU天花板上的高清摄像头捕捉一系列图像以监控植物的状态。为了确保正确的图像曝光，光照强度被减弱（最多2分钟）到显著较低的水平，导致植物的光合作用速率暂时降低。第二个事件是指向营养液充氧，以保持溶解氧浓度高于根系呼吸所需的临界水平。该策略包括从PCU大气循环中向营养液泵送空气，如果（连续测量的）溶解氧浓度降至给定设定点以下。在每个沟槽底部放置的多孔石头中向营养液充入空气，这个过程对于高度发达的根系可能需要几十分钟。在充氧事件期间，观察到PCU大气模块中的CO2浓度短暂增加。这种效应可归因于营养液中的碳酸盐平衡发生偏移，其中溶解CO2的部分压力降低促进了碳酸氢盐向气态CO2的转化。在这种情况下，这种现象由于营养液的亚酸性条件（pH = 5.9）而加剧，导致CO2释放到PCU大气中。因此，这些数据也从用于图8的集合中移除。结果提供了几个值得注意的见解。首先我们注意到，即使使用观察器提供的估计值来增强一个简单的PI也能显著提高整体性能。在第18天之后，跟踪误差不仅大幅减少（减少了三倍），而且对植物生长的敏感性也降低了，这在实验的后期阶段越来越明显。这种改进在图9中可以进一步看到，我们展示了第18天的详细情况，包括切换到更高级实现的时刻。下载：下载高分辨率图像（421KB）下载：下载全尺寸图像图8. 两个生命周期实验中CO2跟踪误差的演变。箱形图的范围从第25百分位数到第75百分位数，它是基于控制器活跃的数据间隔计算的（详见第6节）。在实验1的第18天，ek的计算从(41)切换到(42)。MPC在整个实验2中证明非常有效，无论植物的发育阶段如何，都显示出高可靠性。与实验1相比，MPC显著降低了平均和最大CO2误差，尽管羽衣甘蓝的CO2需求范围（约10%）高于生菜。即使与高级PI策略相比，这种改进也是显著的，平均误差从大约6.6 ppm降低到接近1 ppm。图10更详细地比较了两种情况下的CO2动态，我们展示了实验2第24天的最后4个小时和实验2第21天的情况。在这些日子里，两个试验中的植物在耗尽阶段的吸收率相似，分别为0.0275 mmol s-1和0.0286 mmol s-1。MPC策略的可靠性显而易见。更重要的是，最小误差使得能够清楚地区分影响CO2水平的其他现象，例如在光采集期间光照强度的降低和曝气器的激活。尽管这些效应在实验1中也存在，但由于高级PI控制器引入的误差幅度相当，因此无法区分。下载：下载高分辨率图像（251KB）下载：下载全尺寸图像图9.实验1在第18天的控制性能细节，重点展示了我们从基于（41）的PI实现切换到基于（42）的实现的时刻。尽管1000 ppm在商业和研究应用中是常见的浓度（Peiro等人，2020年；Zhang等人，2023年），但PCU仍在超高的CO2浓度下运行，以重现封闭空间栖息地典型的关键大气情景。这些试验的主要目的是评估植物在极端CO2暴露条件下的表现，这些条件与BLSS应用相关。同时，这些测试还验证了PCU硬件和控制软件在高浓度运行状态下的性能。使用表1中实验2相同的作物和环境配置，进行了两个额外的体内试验，设定点分别为5000 ppm和10000 ppm。图11展示了这三种情况下绝对跟踪误差|κ??κk|的统计分析。结果表明，即使设定点比标准值高一个数量级，调节精度仍然得以保持，达到了测量仪器几乎无法检测到的浓度水平。值得注意的是，在10000 ppm的实验过程中，误差始终低于使用先进PI控制器在1000 ppm时观察到的误差。

图10. 实验1（左）和实验2（右）第24天最后4小时的光照条件下CO2调节细节（黑暗阶段转换在第16小时）。选择这些具体日期是为了确保两个实验之间的CO2吸收量相匹配。底部图表显示了光照配置（1对应于低强度辐射以进行图像采集，0则表示其他情况）和曝气器激活状态。

尽管实验1和实验2涉及的作物不同（表1），但两者都在相同的环境条件下进行，并且都处于生长期。控制问题与具体作物无关：腔室动态仅取决于PCU内部的CO2交换总量，这包括植物的净交换量和其他未建模的贡献（例如，微生物活动或营养液中的CO2释放）。因此，不同作物仅通过这种缓慢变化的干扰的幅度和变异性影响系统，而提出的估计器-控制器框架旨在实时处理这些变化。为了提供一个严格的控制器比较基准，我们报告了另一个实验，在该实验中，控制器在同一天内依次在相同的作物上进行了评估。选择蚕豆（Vicia faba L.）是因为它在BLSS中作为主要蛋白质作物的战略作用以及其对于长期太空任务的高营养价值。蚕豆的显著生物量生产和明显的代谢活动导致CO2吸收需求增加，为评估封闭空间环境中高代谢负荷下的控制器性能提供了代表性的测试案例。在试验期间，MPC和先进的PI控制器在相同的环境条件下依次应用，从而允许进行直接和无偏的性能比较。图12详细展示了实验第42天从MPC控制器到先进PI控制器的转换过程。在这个阶段，植物吸收量估计为0.053 mmol s-1——几乎是实验1和实验2观察到的速率的两倍。这一测试证实了我们之前的发现，即MPC的平均误差为1.32 ppm，而PI的误差为7.69 ppm。与图10的结果类似，较小的调节误差有助于表征影响CO2水平的现象。重要的是，我们强调所有实验中都使用了一组相同的校准参数[w,p]。相比之下，由于植物动态的影响加倍，PI控制器在这次最终实验中需要重新调整以达到与早期测试相当的性能。

图11. 在表1中与实验2相同的作物和环境条件下，三个体内试验中CO2跟踪误差的演变，设定点分别为1000 ppm、5000 ppm和10000 ppm。

图12. 在蚕豆（Vicia faba L.）栽培的体内实验中，从第1天到第42天的跟踪误差。图表突出了所提出的MPC控制器与第3节中增强的PI控制器之间的切换，显示了MPC获得了更准确的调节效果。

7. 讨论与可能的扩展
我们首先强调，高跟踪精度对于PCU至关重要，主要是因为它在收集研究应用（如模型合成）所需的数据方面的作用。如图10所示，当跟踪误差较大时，某些效应的动态难以观察（因此也难以准确建模）。然而，在植物生长需求的背景下，使用先进PI实现所达到的约10 ppm的CO2跟踪误差可能已经足够准确。这表明，本文提出的方法即使在要求较低的情景中也具有价值，因为估计器显著提高了基本控制解决方案的性能，使其在作物变化的情况下也能足够准确。在其他情境中，例如农业或BLSS集成中，控制目标通常是根据操作范围而不是固定设定点来定义的。在这种情况下，标准控制器无法利用这种额外的自由度，仍然需要指定一个虚构的设定点（通常选择为范围的中点）。相比之下，MPC可以轻松适应纯粹由CO2约束驱动的需求，从而提高灵活性并促进系统级优化（Cimini等人，2025年）。通过稍微修改（34），我们得到了公式：(43a)min{cj+1,uj}j=0p?1∑j=0p?1J(δj,ej+1)(43d)s.t.cj+1=cj+ts(uj+??j|k)cmin≤ci+1≤cmaxuj∈{0}∪[umin,umax]，其中cmin,cmax是浓度界限，可能是随时间变化的，并由监督者定义。作为最后的说明，这两种公式都可以很容易地扩展以纳入基于已知次要效应的预测控制动作，这些效应会影响?的动态。例如，可以提前轻松预测光照强度或曝气器激活的可预见变化，并利用这些信息更准确地估计?在p步预测范围内的未来演变。这使得可以设计出更简洁的控制器，当预期由于外部事件导致CO2增加时，该控制器可以避免不必要的CO2注入，同时仍保持浓度在规定的范围内。这在为回收优化的BLSS背景下是一个有价值的特性（Cimini等人，2025年）。

8. 结论
我们提出了一种准确调节生长室内CO2浓度的方法。该方法被证明在多种具有挑战性的情景中都是可靠的，而在这些情景中，现有方法通常表现不佳。这些情景包括：操作范围内缺乏或不准确的作物模型、同时种植多种作物、关于作物发育阶段的不确定性，以及灵活或未定义的收获计划。该方法依赖于两个关键组成部分。第一个是一个静态的、与作物无关的卡尔曼滤波器，它估计腔室内的自由CO2交换量，并选择性地更新完整状态向量的子集以处理由光照引起的植物切换动态。该估计器旨在本质上处理浓度测量的典型延迟。第二个是一个混合MPC，旨在应对众所周知的执行器限制，提供更高的控制精度并满足高级控制要求。尽管在实验PCU环境中得到了验证，但所提出的估计器-控制器框架具有与作物无关的结构和计算效率，使其同样适用于大规模受控环境农业，在这种农业中，异质种植方案和灵活的生产计划很常见。该方法的一个副产品是估计植物的光合作用和呼吸活动，我们计划在未来将其用于增强在更复杂的BLSS框架中的集成，以及用于诊断任务。正如我们在实际实验中展示的那样，其性能远优于其他方法，也包括Pannico等人（2022年）最初提出的方法。因此，该方法现在已成为位于那不勒斯费德里科二世大学空间作物研究实验室的ESA-MELiSSA PCU官方软件发布的一部分。未来的实验测试将验证控制算法在极端环境干扰下的鲁棒性。

CRediT作者贡献声明
Gionata Cimini：撰写——原始草稿、软件、方法论、概念化。
Antonio Pannico：撰写——原始草稿、验证、调查、概念化。
Marco Gatti：撰写——审阅与编辑、数据管理。
Daniele Bernardini：撰写——审阅与编辑、监督。
Stefania De Pascale：撰写——审阅与编辑、项目管理、资金获取。

资金支持
该研究得到了欧洲航天局（ESA）和法国的支持，通过MELiSSA项目PlAnt封闭生命支持系统——工程、制造和测试（PaCMAN）提供资金。