下面我先把你老师提到的“模糊点”逐条对应成一条清晰的技术主线，然后给你一份可以直接替换进论文的**“更清楚的第二～四章写法（含变量定义、厚度/张力位置、参数化合理性、网络设计动机、控制目标/参考/冲突/求解、预测与控制的衔接）”**。你可以按你论文的具体变量维度（比如每个机架状态到底包含哪些量）再做小改。

0. 先把整体逻辑理顺（你现在的核心问题）

你现在的写法里，读者会卡在 6 个地方：

1. Δx、x、u、Δu 的物理含义没落地：到底是“偏差量(deviation)”还是“相邻采样差分(increment)”？厚度、张力具体在哪个分量？
2. 式(2)的“区间内多项式参数化”没有解释为什么在五机架冷连轧可行：控制是数字采样的，为什么要搞连续时间 τ？
3. 网络结构没讲清楚“有 u”和“没 u”的差别：MPC 里必须能“随着你改变控制而改变预测”，否则控制优化没有意义。
4. 训练目标的动机不清楚：多步 rollout 为啥要做？forward/backward reciprocal consistency 是解决什么问题？
5. 第四章读者看不出控制策略怎么落地：参考信号是什么、优化目标怎么写、耦合冲突在哪里、怎么协调、最终怎么解、怎么滚动执行。
6. 预测与控制怎么衔接没讲：NN 输出是什么、怎么进入 MPC 预测、怎么在每次采样时闭环更新。

你需要做的不是“再写更多公式”，而是把每一章都改成：先把对象定义清楚 → 再说明为什么这么建模/这么参数化 → 再给算法流程（输入/输出/求解器）。

1) 建议你统一符号：把“Δ”讲清楚（这是最大雷点）

你现在把 Δx 叫“state increment”，但又用它做状态演化。为了让读者不迷路，强烈建议你在论文里二选一并固定：

方案 A（更推荐，最稳）：Δ 表示“偏差量 deviation”

• 令 $x_i(t)$ 是第 i 机架相关物理量（厚度、张力等）的状态向量
• 选一个工作点/名义轨迹 $x_i^{\star}(t)$，定义

同理

这样你的离散映射“Δx → Δx”就完全合理：它是偏差系统。

方案 B（不推荐但可用）：Δ 表示“采样差分 increment”

那就必须写清楚

但这样进入 MPC 会更绕，而且你的残差网络“shortcut 直接加 Δx(t_n)”会更难解释其物理意义。

建议你用方案 A，把“increment”改成“deviation”或“deviation state”，你老师提的“变量含义不清”会立刻好一半。

2) 厚度和张力到底在哪：给一个“机架 i 的状态向量”示例（你必须明确写）

五机架冷连轧典型关注两类输出：

• 出口厚度：每架出口厚度 $h_i$（或其偏差 $\Delta h_i$）
• 机架间张力：第 i 与 i+1 架之间张力 $T_i$（或 $\Delta T_i$）

一种很常见、读者也最好理解的状态选取（仅示例，你按你仿真模型调整维度）：

• 对机架 $i$（1～5），定义局部状态（偏差）

其中

• $\Delta h_i$：第 i 架出口厚度偏差
• $\Delta T_{i-1}$：(i-1,i) 间张力偏差（i=1 时没有，可置 0 或不纳入）
• $\Delta T_i$：(i,i+1) 间张力偏差（i=5 时没有）

邻居集合也要落地写清楚（五机架串联系统最简单）：

3) 式(2)的“区间内二次多项式参数化”在五机架场景为什么合理？（你要这样解释）

你的参数化不是为了“数学好看”，而是为了解决冷连轧控制里两个现实点：

1. 执行器/设定值需要平滑：辊缝/速度/弯辊等设定值如果采样点间跳变，会激发张力振荡、厚度波动。
2. 采样周期内控制并非真正常值：工业上常见做法是采样更新“目标”，但底层驱动/液压闭环会在采样周期内连续变化；用低阶多项式近似一个采样区间内的变化是合理的“等效输入模型”。

因此你可以把区间内的控制偏差写成基函数展开（你这里选二次）：

这等价于“用 3 个参数描述一个采样区间内的平滑控制轨迹”，比 piecewise-constant 更能表达渐变/缓变，也更利于约束“采样区间内任何时刻不越界”。

然后你定义的离散等效输入用“区间平均”：

这一步要强调一句：

控制优化决策变量是 $\Gamma_{i,n}$，但为了与离散状态更新与执行器指令一致，我们用区间平均将连续轨迹映射为离散等效控制量，同时区间内约束仍用轨迹检查保证。

再补一句你老师关心的：“可不可以在五机架情境这样参数化？”
答案是可以，但要加一句适用条件：采样周期相对执行器带宽不能太大、且实际指令在区间内变化可被低阶多项式近似。你可以这样写：

当采样周期 $\delta_n$ 相对液压/速度内环带宽不至于过大时，区间内的控制偏差轨迹可用低阶多项式近似；二次项在不显著增加决策维度的前提下提供“曲率自由度”，能更好拟合加减速/辊缝渐变等工况。

4) “有 u”和“没 u”到底差在哪？（第三章必须讲透）

这句话你可以直接写进论文：

• 如果网络输入里没有控制量（或其参数 $\Gamma$），网络学到的是“在训练数据控制策略下系统如何演化”的闭环轨迹模型。这模型不能用于 MPC，因为 MPC 要通过改变控制来改变预测；没有控制输入的模型无法回答“如果我换一个控制会怎样”。
• 如果网络输入包含 $\Gamma_{i,n}$（以及 $\delta_n$），网络才是“可控预测模型（surrogate dynamics）”，能用于在线优化。

所以你现在的写法必须明确：你训练的是

而不是“仅靠状态自回归”。

5) 为什么用 residual 结构？为什么加 η 分支？（你要给明确动机）

你现在写得像“为了鲁棒性加个分支”，但读者不知道它解决啥。建议这样讲：

residual（shortcut）解决什么

• 轧机采样周期一般不大，相邻采样时刻状态偏差变化往往“相对平缓”。
• 用残差结构相当于给了一个物理先验：短时内状态变化不大（persistence prior）。
• 网络只学习“校正项”，训练更稳定、多步滚动误差更小。

变采样/大 δ 的问题

当 $\delta_n$ 变大时，系统变化幅度变大，单一残差网络容易对“时间尺度”敏感。你引入

一个更清晰的解释是：

• $\eta_i$ 用来学习随 $\delta_n$ 缩放/偏置的低频项（例如近似线性随 $\delta$ 放大），
• $\mathcal I_i$ 学习剩余的非线性耦合细节。
  这样做能缓解“不同采样间隔下同一状态变化尺度不同”的问题。

6) 第四章你要把“控制策略怎么落地”写成可执行算法

你老师说“没看出来控制策略怎么设计”，通常是因为你缺了 5 件事：

1. 控制目标到底是什么：厚度跟踪？张力稳定？两者权衡？
2. 参考信号是什么：常值设定？来自工艺计划？由上层调度给？
3. 耦合冲突在哪里：一个机架改速度会影响两侧张力，厚度与张力可能相互拉扯。
4. 协调机制是什么：你说 Nash，但要说清楚每个子问题里“邻居的信息”怎么进入预测和代价。
5. 最终怎么解：每个子问题是 NLP 还是 QP？用什么求解器？每次迭代怎么停止？滚动时怎么应用？

下面给你一份“更清楚的写法模板”。

✅ 给你一份可直接替换论文的“重写版思路 + LaTeX骨架”

你可以把下面当成第二～四章的“清晰版本”，再按你实际状态/控制维度微调。

第二章（建议重命名）数据集构造与离散等效输入

2.1 五机架耦合系统与符号定义（必须新增）

（建议放一个符号表/文字说明）

• 五机架编号 $i\in\{1,2,3,4,5\}$，采样时刻 $t_n$，区间长度 $\delta_n=t_{n+1}-t_n$。
• 机架 $i$ 的局部状态选为（示例）

其中 $h_i$ 为第 i 架出口厚度，$T_i$ 为 (i,i+1) 间张力；$(\cdot)^\star$ 为名义工况/参考轨迹。

• 控制输入（示例）可以是辊缝设定、速度设定或其组合：

若只选单输入（例如辊缝），则 $u_i\in\mathbb R$。

• 邻居集合（串联）：

并定义邻居堆叠

这一段写清楚后，厚度和张力“在哪”就明确了：它们就是 $x_i$ 的分量。

2.2 离散映射与等效输入的定义（把你式子“物理化”）

由于五机架系统存在显著耦合、摩擦与材料参数不确定等因素，难以获得高精度的一阶机理离散模型。本文采用数据驱动方式学习如下局部离散演化关系：

其中 $\Phi_i(\cdot)$ 为未知非线性映射；$\Delta u_i([t_n,t_{n+1}])$ 表示该采样区间内的控制轨迹信息。

为与离散预测/执行一致，引入区间平均等效控制偏差：

相应地，扰动也可用区间平均表示 $\Delta d_i(t_n)$。

2.3 区间内控制轨迹参数化（回答“能不能这么做”）

为在不显著增加在线优化维度的前提下描述采样区间内的平滑控制变化，本文将区间控制偏差轨迹用二次多项式参数化：

其中 $\Gamma_{i,n}=[\gamma_{i,n0},\gamma_{i,n1},\gamma_{i,n2}]^\top$。
该参数化可近似描述轧机执行器在采样周期内的渐变响应，并天然保证控制轨迹在区间内可检查约束（见第4章）。

由此区间平均等效输入为

（这里加一句适用条件/解释：采样周期与执行器带宽关系、平滑性需求。）

2.4 一步样本与 K 步段样本（把“输入/输出”写成机器学习可读格式）

一步样本（监督学习）
对每个区间 $[t_n,t_{n+1}]$，构造训练输入

输出标签为 $\Delta x_i(t_{n+1})$ 或残差 $\Delta r_i=\Delta x_i(t_{n+1})-\Delta x_i(t_n)$。

K 步段样本（用于 rollout 与一致性正则）
组织连续 K 个区间得到段样本 $\mathcal W_{i,n}$，包含

这一段的意义要写明：用于训练时约束长期递推误差（drift），让模型可用于 MPC 的长预测域。

第三章 残差神经网络预测模型（必须把“有u/无u、为什么残差”讲透）

3.1 预测任务定义：我们要学什么？

目标是学习一个可控的一步预测模型：

强调：$\Gamma_{i,n}$ 进入输入，保证模型能响应控制变化，从而可用于 MPC 在线优化。

3.2 残差结构与网络设计动机

• shortcut 项 $\Delta x_i(t_n)$ 提供“短时保持”的先验；
• 网络仅学习校正项（包含耦合与非线性），多步滚动更稳定；
• 引入 $\delta_n$ 使模型具备变采样泛化能力。

你原来的 $\eta_i + \mathcal I_i$ 可以这样解释并写清楚：

其中 $\eta_i$ 用于学习与时间尺度 $\delta_n$ 强相关的低频项（如随 $\delta$ 的尺度变化），$\mathcal I_i$ 学习剩余耦合非线性细节；当 $\eta_i\equiv 0$ 退化为标准 ResNet。

3.3 多步训练：rollout loss + reciprocal consistency 为何需要？

要用于 MPC，模型必须在 $N_p$ 步滚动中不发散。仅用一步 MSE 容易出现“短期准、长期漂”。
因此引入：

• rollout loss：直接惩罚多步预测轨迹与真值偏差（抑制误差累积）
• forward/backward reciprocal consistency：通过学习一个反向模型 $\mathcal B_i$ 并约束“正向滚到终点再反推应回到起点”，从结构上约束模型减少漂移（类似可逆一致性的正则化思想）

（这段你现在写了很多公式，但缺“这一招解决什么问题”的解释。补上这两三句话，老师就容易认可。）

第四章 RNE-DMPC 控制设计（你要把“参考、冲突、求解、滚动执行”写成完整闭环）

4.1 控制目标、参考信号与冲突来源（必须写清楚）

• 控制目标：

  1. 厚度跟踪：使各机架出口厚度偏差 $\Delta h_i$ 跟踪参考 $\Delta h_i^{\text{ref}}$（通常为 0 或来自工艺设定）
  2. 张力调节：使机架间张力偏差 $\Delta T_i$ 跟踪 $\Delta T_i^{\text{ref}}$（通常为 0 或张力设定偏差）

• 参考信号来源（你必须选一种说法并固定）：

  • 常值设定：$\Delta x_{\text{ref}}(t)\equiv 0$（围绕名义工况稳态调节）
  • 或上层工艺给定的时变轨迹：$\Delta x_{\text{ref}}(t_{n+s})$ 来自生产计划/换规格策略/张力与厚度设定曲线

• 冲突来源（耦合）：
  在五机架串联系统中，机架 $i$ 的速度/辊缝等动作会改变带材流动与变形，从而影响相邻机架间张力；即 “局部最优”可能使邻居张力恶化。
  这就是为什么需要分布式协调（Nash 迭代）而不是各做各的。

4.2 用 NN 做预测：把“预测如何服务控制”写成一句话

在 MPC 每次在线优化中，候选决策变量是 $\mathbf\Gamma_i$。给定 $\mathbf\Gamma_i$ 和邻居的最新策略/预测，利用已训练网络递推得到预测轨迹：

这条递推就是“控制优化的模型约束”。 没这条，MPC 就无从谈起。

4.3 局部优化问题：决策变量、代价函数、约束（写清楚“是什么 + 为什么”）

• 决策变量（控制域 $N_c$）：

• 代价函数（把厚度/张力写进来，别只写 $\Delta x$）：

若 $\Delta x_i=[\Delta h_i,\Delta T_{i-1},\Delta T_i]^\top$，则可写

（你可以把“厚度权重更大/张力权重更大”的工艺解释写一句。）

• 约束（要说清楚“如何把 Γ 变成对 u 的约束”）：

区间内控制轨迹约束：

实现上检查 $\tau=0,\delta$ 以及二次极值点 $\tau^\star=-\gamma_1/(2\gamma_2)$（若落在区间内）。

绝对输入约束通过区间平均更新：

从而 enforce：

4.4 Nash 协调：冲突怎么解决？（用“最佳响应”把逻辑讲通）

• 每个子系统 $i$ 在给定邻居策略 $\mathbf\Gamma_{Z_i}$ 下求解自己的最优响应（best response）：

• 通过通信交换 $\mathbf\Gamma_i^{(l)}$ 与预测轨迹 $\Delta\hat x_i^{(l)}$，迭代直到收敛：

这段要明确说：Nash 迭代的本质是**“耦合冲突通过重复互相回应找到一个稳定折中点”**，它不是强行让大家一致（那是共识/ADMM），而是让谁也不愿意单方面改变（均衡）。

4.5 最终怎么解？（你必须给出“求解器层面”的一句话）

每个局部问题通常是一个带 NN 预测约束的非线性规划（NLP）。因为 NN 可微，可以采用：

• SQP / interior-point（如 IPOPT）
• 或基于自动微分的梯度法（若你用 PyTorch/TF 端到端）

你不需要写“我用 IPOPT”，但至少要写：

本文将局部优化问题视为可微 NLP，利用梯度信息进行数值求解；每次 Nash 迭代求得 $\mathbf\Gamma_i^{(l)}$ 后广播更新。

4.6 滚动时域闭环：控制与预测如何衔接（用 5 行把闭环讲清）

在每个采样时刻 $t_n$：

1. 测量/估计当前 $\Delta x_i(t_n)$，并接收邻居上一次迭代/上一时刻的预测或策略初值；
2. 运行 Nash 迭代：每个机架并行求解局部 NLP，交换预测与策略，直至收敛；
3. 仅实施第一段控制参数 $\Gamma_{i,n}^\star$，生成区间内控制轨迹 $\Delta u_{i,n}(\tau)$；
4. 用区间平均更新实际输入 $u_i(t_n)$，系统运行到 $t_{n+1}$；
5. $n\leftarrow n+1$ 滚动重复。

这段写出来，读者就能看出“预测—优化—执行—再测量”的闭环链条，你老师说的“控制和预测怎么衔接”就解决了。

7) 你可以直接照着改的“最小清单”（不加新内容也能显著变清楚）

如果你不想大改结构，至少做这 10 个动作：

1. 在第二章开头新增“符号与状态/输入定义”（写清厚度、张力在哪个分量）。
2. 把“increment”改成“deviation”或明确 Δ 的定义（二选一）。
3. 在式(2)后加一段解释：为什么区间参数化、为什么二次、适用条件。
4. 在数据集定义里明确：网络输入是 $[\Delta x_i,\Delta x_{Z_i},\Gamma,\delta]$，输出是 $\Delta x_i(t_{n+1})$ 或残差。
5. 第三章新增一段：“没有 u 的模型不能做 MPC，有 u 的模型才可控预测”。
6. 第三章新增一段：残差结构为何能稳定多步、为何引入 η 分支应对 δ 变化。
7. 第三章训练目标每个 loss 前加一句“它解决什么”。
8. 第四章开头先写“控制目标=厚度+张力；参考=0 或工艺给定”。
9. 第四章明确“冲突=耦合张力传播导致局部动作影响邻居”。
10. 第四章给出“在线闭环流程（5步）+ 每个子问题是 NLP，用可微优化求解”。

如果你愿意，我也可以在你“真实的状态向量/控制量定义”（比如你到底选了辊缝、速度、弯辊还是张力给定等）基础上，把上面的“示例版”完全改成与你论文一致的最终可直接粘贴版：你只要告诉我两件事就行——

1. 你每个机架的 $x_i$ 具体包含哪些量（厚度、张力之外还有没有速度差、压下量、轧制力等）以及维度 $d$；
2. 你的 $u_i$ 到底是哪一个/哪些执行量（单输入还是多输入）。

在你这个“五机架冷连轧、状态只含厚度+张力”的设定下，我对 $u_i$ 的理解是：

$u_i$ 是机架 $i$ 在在线控制中“能直接下发/能被控制器改变的执行量（manipulated variables）”的集合，也就是你 MPC 的决策变量最终要落到的“设定值/指令”。

在冷连轧里，能同时影响厚度与张力、且工业上最常用、最容易让审稿人/老师接受的选择通常是下面这两个：

1) 最推荐、最标准的解释：$u_i=[s_i,\; v_i]$

(1) 辊缝/压下（Hydraulic gap / screw-down）

• 记作 $s_i$（辊缝设定或等效压下量，很多文献也用 $S_i$、$h_{0,i}$、$x_{g,i}$）
• 主要影响：机架出口厚度 $h_i$（“厚度控制”最直接的操纵量）
• 物理直觉：辊缝变小 → 压下率增大 → 出口厚度变薄（当然具体关系受轧制力、材料硬化等影响，所以你用数据驱动去学）

(2) 机架速度（roll speed / strip speed reference）

• 记作 $v_i$（电机速度/线速度设定）
• 主要影响：张力（张力来源于相邻机架速度差 + 带钢延伸/弹性效应）
• 物理直觉：若 $v_{i+1}$ 相对 $v_i$ 增大，(i,i+1) 间张力 $T_i$ 往往上升（方向取决于你张力定义与模型）

所以写成：

这样你第四章要做“厚度–张力协调控制”，逻辑就非常顺：

• 厚度靠 $s_i$ 纠偏
• 张力靠 $v_i$ 协调
• 耦合冲突来自“速度改动会影响上下游张力与厚度演化、辊缝改动也会影响带钢流动从而间接影响张力”。

这也能直接回应老师质疑：“有 $u$ 和没 $u$ 的区别”：没 $u$ 你就根本没法通过优化改变预测张力/厚度。

2) 你现在的“区间二次多项式参数化”如何对应到多输入？

如果 $u_i$ 是二维（辊缝+速度），那你现在的参数化只要“对每个输入通道各来一套二次系数”即可：

于是

• 决策向量维度从你原来“单输入 $p=3$”变成“双输入 $p=6$”：

• 区间平均也分别计算（对每个通道同样是 $+\delta/2,+\delta^2/3$ 那套）

这点你在论文里一定要一句话点明：

“当 $n_u>1$ 时，$\Gamma_{i,n}$ 由各输入通道的多项式系数组合而成，维度为 $p=3n_u$。”

3) 另一种也能自洽、但需要你说明“内环”的选择：$u_i=[s_i,\; T_i^{ref}]$

有些产线确实存在“张力内环”（通过调速实现张力跟踪），这时外环 MPC可以不直接给 $v_i$，而是给一个张力参考 $T_i^{ref}$：

• $u_i=[s_i,\;T_i^{ref}]$
• 速度由底层张力控制器自动调整

但注意：如果你这么定义 $u_i$，你的数据集离线仿真必须把“张力内环控制器”也包含进被控对象，否则网络学到的映射和现场执行不一致。老师如果问“你到底控制的是速度还是张力参考”，你必须说清楚层级结构。

4) 你现在只写“厚度+张力”，那 $u_i$ 选什么更合理？

在“只关心厚度+张力”且你要做 RNE-DMPC（分布式协调） 的语境里，我建议你论文里直接采用最经典的一句：

每机架采用两个可操纵量：辊缝（压下）与速度设定。辊缝主导厚度，速度差主导机架间张力；两者通过带钢变形与张力传播产生耦合，因此需要分布式协调。

这样第二章的参数化、第三章的 NN 输入、第四章的优化目标/冲突来源都能一条线串起来。

如果你愿意，我可以按你“状态只含厚度+张力”的约束，帮你把文中所有出现的 $u_i$、$\Delta u_i$、$\Gamma$、$p$ 以及第4章代价函数里“$\Delta x_{\rm ref}$”改写成明确的厚度参考 $\Delta h^{ref}$ 和 张力参考 $\Delta T^{ref}$，并把“冲突是什么、Nash 怎么解决、每个子问题怎么解、滚动怎么执行”那一段写成你论文可直接替换的 LaTeX 段落。你只要告诉我你希望：

• 每个机架的状态 $\Delta x_i$ 具体取 $[\Delta h_i,\Delta T_{i-1},\Delta T_i]$ 还是更简（例如只取 $[\Delta h_i,\Delta T_i]$）；
• 张力是定义为机架间 $T_i$（i=1..4）还是每机架一个张力量。

latex
%========================
\section{Construction of Dataset}
%========================

\subsection{Five-stand tandem mill setting and notation (variables, thickness and tension)}
\label{subsec:notation_5stand}

Consider a five-stand tandem cold rolling mill indexed by $i\in\{1,2,3,4,5\}$.
Let $t_n$ denote the $n$-th sampling instant and $\delta_n=t_{n+1}-t_n$ the $n$-th sampling interval length.

\textbf{Thickness and tension variables.}
Let $h_i(t)$ denote the exit thickness of stand $i$, and let $T_i(t)$ denote the inter-stand tension between stand $i$ and $i+1$
(i.e., $T_i$ is defined for $i=1,2,3,4$). The control goal of this paper is to achieve coordinated regulation/tracking of
$\{h_i\}_{i=1}^5$ and $\{T_i\}_{i=1}^4$ in the presence of strong inter-stand coupling.

\textbf{State (deviation) definition.}
To avoid ambiguity, throughout this paper the symbol ``$\Delta$'' attached to a \emph{state} denotes a \emph{deviation} (tracking error)
with respect to a prescribed reference (or nominal) trajectory:
\begin{equation}
\Delta h_i(t)\triangleq h_i(t)-h_i^{\mathrm{ref}}(t),\qquad
\Delta T_i(t)\triangleq T_i(t)-T_i^{\mathrm{ref}}(t).
\end{equation}
For the five-stand case, we choose a stand-wise local state vector that explicitly contains the thickness and the adjacent tensions:
\begin{equation}
\Delta x_i(t)\triangleq
\begin{bmatrix}
\Delta h_i(t)\\
\Delta T_{i-1}(t)\\
\Delta T_i(t)
\end{bmatrix}\in\mathbb{R}^{d},\qquad d=3,
\label{eq:local_state_def}
\end{equation}
with the boundary convention $\Delta T_{0}(t)\equiv 0$ and $\Delta T_{5}(t)\equiv 0$ so that all stands share a unified dimension $d=3$.
In the cost function and constraints, boundary ``virtual'' tensions can be assigned zero weights so that they do not affect the optimization.

\textbf{Neighbor set (five-stand chain coupling).}
The mill coupling is dominated by tension propagation between adjacent stands, therefore we use the chain neighbor sets
\begin{equation}
Z_1=\{2\},\quad
Z_i=\{i-1,i+1\}\ (i=2,3,4),\quad
Z_5=\{4\}.
\label{eq:neighbor_set}
\end{equation}
Define the neighbor-state stack (used to encode inter-stand coupling information) as
\begin{equation}
\Delta x_{Z_i}(t_n)=\mathrm{col}\{\Delta x_k(t_n)\,|\,k\in Z_i\}.
\label{eq:neighbor_stack}
\end{equation}

\textbf{Control inputs (actuators) and control increments.}
For each stand $i$, the manipulated variables are chosen as the \emph{roll gap (screw-down/hydraulic gap)} and the \emph{stand speed}:
\begin{equation}
u_i(t)=
\begin{bmatrix}
s_i(t)\\
v_i(t)
\end{bmatrix}\in\mathbb{R}^{n_u},\qquad n_u=2.
\label{eq:ui_def}
\end{equation}
To ensure smooth actuator evolution and to match typical industrial implementations, the optimization is conducted on
\emph{control increments} (sample-to-sample changes) rather than absolute inputs:
\begin{equation}
\Delta u_i(t_n)\triangleq u_i(t_n)-u_i(t_{n-1})
=
\begin{bmatrix}
\Delta s_i(t_n)\\
\Delta v_i(t_n)
\end{bmatrix}.
\label{eq:du_increment_def}
\end{equation}
Note the distinction: $\Delta x$ is a \emph{deviation state} (tracking error), while $\Delta u$ is an \emph{input increment}.

\textbf{Disturbance.}
Let $d_i(t)$ denote exogenous disturbances (e.g., entry thickness fluctuation, friction variation, material property changes, etc.).
We use $\Delta d_i(t)$ to denote their interval-level equivalent representation (defined later).

\subsection{Discrete interval mapping and data-driven motivation}
\label{subsec:discrete_mapping}

The stand-wise deviation-state evolution over $[t_n,t_{n+1}]$ can be expressed by an equivalent discrete-time mapping
\begin{equation}
\Delta x_i(t_{n+1})
=
\Phi_i\Big(\Delta x_i(t_n),\,\Delta x_{Z_i}(t_n),\,\Delta u_i([t_n,t_{n+1}]),\,\Delta d_i([t_n,t_{n+1}])\Big),
\label{eq:true_unknown_mapping}
\end{equation}
where $\Phi_i(\cdot)$ is generally nonlinear and strongly coupled due to tension propagation and rolling deformation interactions.
A commonly used \emph{conceptual} local linear discrete form is
\begin{equation*}
\Delta x_i(t_{n+1})
=
M_d\,\Delta x_i(t_n)
+
N_d\,\Delta u_i(t_n)
+
F_d\,\Delta d_i(t_n),
\end{equation*}
where $M_d,N_d,F_d$ are equivalent discrete matrices and $\Delta u_i(t_n),\Delta d_i(t_n)$ are interval-level equivalent inputs/disturbances.
However, in a practical five-stand cold rolling mill, accurate identification of $(M_d,N_d,F_d)$ from first principles is difficult,
because of complex coupling, unmodeled nonlinearities, and time-varying operating conditions.

\begin{remark}
In fact, due to the existence of complex coupling relationships, it is difficult to directly and accurately establish the discrete mapping
based on first principles. Therefore, in this paper, we aim to learn an accurate approximation of \eqref{eq:true_unknown_mapping} from data,
and then embed the learned surrogate into distributed MPC.
\end{remark}

%========================
\subsection{Interval-level parameterization and one-step dataset}
%========================

To construct the supervised dataset for training, we locally parameterize the \emph{control increment trajectory} within each sampling interval.
For the local interval $[t_n,t_{n+1}]$, define the interval length
\begin{equation}
\delta_n = t_{n+1}-t_n ,
\end{equation}
and introduce a local time variable $\tau\in[0,\delta_n]$.

\textbf{Why interval-level parameterization is valid in five-stand cold rolling.}
Although the controller updates the setpoints at discrete instants $t_n$, the physical actuators (hydraulic gap and drive speed loops)
evolve continuously within $[t_n,t_{n+1}]$ and are typically implemented with smoothing/ramps.
Moreover, abrupt changes of roll gap/speed can excite tension oscillations and degrade thickness stability.
Therefore, representing the within-interval increment trajectory by a low-order polynomial provides:
(i) a compact finite-dimensional decision variable for optimization;
(ii) a smooth within-interval command profile;
(iii) a convenient way to enforce \emph{continuous-time} bounds in the whole interval.
This approximation is appropriate when the sampling interval is not excessively large compared to actuator bandwidth and
the within-interval evolution can be well captured by a low-order basis.

\textbf{Vector quadratic polynomial parameterization (two inputs).}
Using a quadratic polynomial basis for each input channel, the control increment trajectory on $[t_n,t_{n+1}]$ is parameterized as
\begin{equation}
\Delta u_{i,n}(\tau;\Gamma_{i,n})
=
\Gamma_{i,n0}+\Gamma_{i,n1}\tau+\Gamma_{i,n2}\tau^2,
\qquad \tau\in[0,\delta_n],
\label{eq:du_poly_vector}
\end{equation}
where $\Gamma_{i,n0},\Gamma_{i,n1},\Gamma_{i,n2}\in\mathbb{R}^{n_u}$ are vector coefficients.
Equivalently, for $n_u=2$ (roll gap and speed), one may write component-wise:
\begin{equation}
\begin{aligned}
\Delta s_{i,n}(\tau) &= \gamma^{(s)}_{i,n0}+\gamma^{(s)}_{i,n1}\tau+\gamma^{(s)}_{i,n2}\tau^2,\\
\Delta v_{i,n}(\tau) &= \gamma^{(v)}_{i,n0}+\gamma^{(v)}_{i,n1}\tau+\gamma^{(v)}_{i,n2}\tau^2.
\end{aligned}
\end{equation}
Define the local parameter vector by stacking all coefficients:
\begin{equation}
\Gamma_{i,n}\triangleq
\big[
(\Gamma_{i,n0})^\top,\,
(\Gamma_{i,n1})^\top,\,
(\Gamma_{i,n2})^\top
\big]^\top
\in\mathbb{R}^{p},
\qquad
p=3n_u=6.
\label{eq:Gamma_dim}
\end{equation}
Here, $\Gamma_{i,n0}$ denotes the initial baseline of the increments, while $\Gamma_{i,n1}$ and $\Gamma_{i,n2}$ describe the linear and quadratic variation rates.

\textbf{Equivalent discrete-time (interval-averaged) increments.}
Define the equivalent discrete-time increments as the interval averages:
\begin{equation}
\begin{aligned}
\Delta u_i(t_n) &\triangleq \frac{1}{\delta_n}\int_0^{\delta_n}\Delta u_{i,n}(\tau)\,d\tau,\\
\Delta d_i(t_n) &\triangleq \frac{1}{\delta_n}\int_0^{\delta_n}\Delta d_i(\tau)\,d\tau.
\end{aligned}
\label{eq:interval_average_def}
\end{equation}
With \eqref{eq:du_poly_vector}, the interval average admits a closed form:
\begin{equation}
\Delta u_i(t_n)=
\Gamma_{i,n0}
+\Gamma_{i,n1}\frac{\delta_n}{2}
+\Gamma_{i,n2}\frac{\delta_n^2}{3}.
\label{eq:interval_average_closedform}
\end{equation}

\textbf{One-step sample generation (five-stand coupled simulation).}
Given the above parameterization, one training sample is generated on each interval $[t_n,t_{n+1}]$.
In addition to the local state deviation, the neighbor state deviations are also included to represent inter-stand coupling.
The specific process is shown in Table~\ref{tab:interval_sample_generation_en}.

\begin{table}[t]
\centering
\small
\renewcommand{\arraystretch}{1.15}
\caption{Procedure for generating one interval-level sample on $[t_n,t_{n+1}]$ (five-stand coupled mill).}
\label{tab:interval_sample_generation_en}
\begin{tabularx}{\linewidth}{>{\centering\arraybackslash}p{0.09\linewidth} X}
\toprule
\textbf{Step} & \textbf{Operation} \\
\midrule
1 & \textbf{State sampling:} sample $\Delta x_i(t_n)$ and its neighbor stack $\Delta x_{Z_i}(t_n)$ from $\mathcal{I}_x$. \\
2 & \textbf{Parameter sampling:} draw $\Gamma_{i,n}\sim\mathcal{I}_\Gamma$ (ranges of polynomial coefficients for both $\Delta s_{i,n}(\tau)$ and $\Delta v_{i,n}(\tau)$). \\
3 & \textbf{Control construction:} compute $\Delta u_{i,n}(\tau)$ via the vector polynomial model \eqref{eq:du_poly_vector}. \\
4 & \textbf{State propagation:} integrate the \emph{five-stand coupled} mill model on $[t_n,t_{n+1}]$ (e.g., RK4) using the within-interval control trajectory, and record $\Delta x_i(t_{n+1})$. \\
\bottomrule
\end{tabularx}
\end{table}

Accordingly, an interval sample for subsystem $i$ can be represented as
\begin{equation}
\mathcal{D}_{i,n}=\big\{\Delta x_i(t_n),\ \Delta x_{Z_i}(t_n),\ \Delta u_{i,n}(\tau),\ \Delta x_i(t_{n+1})\big\},
\label{eq:interval_sample_def}
\end{equation}
which is used to learn the mapping relationship from the current local and neighbor deviation states and the local control-increment trajectory
to the next local deviation state.

For each subsystem $i$, by repeating the above procedure across multiple intervals and randomized draws,
the local one-step training dataset is formed as
\begin{equation}
\begin{split}
S_i=\Big\{&
\big(\Delta x_i^{(j)}(t_n),\,\Delta x_{Z_i}^{(j)}(t_n),\,\Delta x_i^{(j)}(t_{n+1});\\
&\qquad \Gamma_{i,n}^{(j)},\,\delta_n^{(j)}\big)
\ \Big|\ j=1,\ldots,J
\Big\}.
\end{split}
\label{eq:one_step_dataset}
\end{equation}
The overall dataset for the five-stand mill can be denoted as $\{S_i\}_{i=1}^{5}$.
The point-cloud visualization of the training dataset is shown in Figure~\ref{2}.

\begin{figure*}[htbp]
  \centering
  \includegraphics[scale=0.5]{picture/Fig2.pdf}
  \caption{Point cloud map of the training dataset.}\label{2}
\end{figure*}

%========================
\subsection{Multi-step rollout segment dataset}
%========================

The one-step set $S_i$ is sufficient for one-step regression, but it is not sufficient for training with multi-step rollout loss
and reciprocal-consistency regularization, because these objectives require ground-truth state trajectories over a horizon of $K$ consecutive intervals.
Therefore, without changing the single-interval sampling mechanism above, we additionally organize the offline-simulated samples
into $K$-step trajectory segments.

Specifically, during offline simulation, for each starting time $t_n$ we generate a segment of length $K$ by consecutively sampling
$\{\Gamma_{i,n+s},\delta_{n+s}\}_{s=0}^{K-1}$ (and the corresponding inputs/disturbances),
and integrating the five-stand coupled mill model over $[t_{n+s},t_{n+s+1}]$ for $s=0,\ldots,K-1$.
Hence, we obtain the deviation-state sequence $\{\Delta x_i(t_{n+s})\}_{s=0}^{K}$ as well as the neighbor stacks
$\{\Delta x_{Z_i}(t_{n+s})\}_{s=0}^{K}$.

Define a $K$-step segment sample for subsystem $i$ as
\begin{equation}
\begin{aligned}
\mathcal{W}_{i,n}=
\Big\{&
\big(\Delta x_i(t_{n+s}),\,\Delta x_{Z_i}(t_{n+s}),\,\Gamma_{i,n+s},\,\delta_{n+s}\big)_{s=0}^{K-1}; \\
&\big(\Delta x_i(t_{n+s+1})\big)_{s=0}^{K-1}
\Big\}.
\end{aligned}
\label{eq:Kstep_segment_def}
\end{equation}
By repeating the above segment generation, we form the multi-step training set
\begin{equation}
S_i^{(K)}=\Big\{\mathcal{W}_{i,n}^{(j)}\ \Big|\ j=1,\ldots,J_K\Big\}.
\label{eq:Kstep_dataset}
\end{equation}
Note that $S_i$ can be viewed as the marginal one-step projection of $S_i^{(K)}$ (by keeping only $s=0$),
thus the original dataset design is preserved, and only an additional \emph{segment organization} is introduced for multi-step training.

%========================
\section{Construction of Residual Neural Network}
%========================

\subsection{Network Architecture (what is learned, why residual, why include $u$)}
\label{subsec:net_arch}

Given the training dataset, the neural network model is defined and trained.
The network model aims to learn the stand-wise deviation-state evolution of the interconnected five-stand cold rolling system over each local interval.
Specifically, for subsystem $i$, the proposed residual network defines a nonlinear mapping
\begin{equation}
\mathcal{N}_i:\mathbb{R}^{d(1+|Z_i|)+p+1}\rightarrow\mathbb{R}^{d},
\label{eq:Ni_map}
\end{equation}
where $d$ denotes the dimension of the local deviation state $\Delta x_i$ ($d=3$ in \eqref{eq:local_state_def}),
$|Z_i|$ is the number of neighbors of subsystem $i$ defined in \eqref{eq:neighbor_set},
and $p=\mathrm{dim}(\Gamma_{i,n})$ denotes the dimension of the local input-parameter vector ($p=6$ for two-input quadratic parameterization).

\textbf{Why the input must include control information.}
If the network input does \emph{not} include the control variables (here, $\Gamma_{i,n}$ and $\delta_n$),
the learned model degenerates to a purely autoregressive predictor that only reproduces trajectories under the \emph{training policy}
and cannot answer ``what will happen if we change the control?''
In MPC, the optimizer must evaluate the predicted trajectory under \emph{candidate} decision variables,
therefore a \emph{control-dependent} predictor (including $\Gamma_{i,n}$) is necessary for online optimization.

Accordingly, we define the residual mapping
\begin{equation}
\begin{aligned}
\Delta r_i &= \mathcal{N}_i(X_{i,\text{in}}; \Theta_i),\\
X_{i,\text{in}} &\in \mathbb{R}^{d(1+|Z_i|)+p+1},\qquad
\Delta r_i \in \mathbb{R}^d ,
\end{aligned}
\label{eq:residual_def}
\end{equation}
where $\Theta_i$ represents the trainable parameters of the network for subsystem $i$, and $X_{i,\text{in}}$ denotes the concatenated input vector for subsystem $i$ in the current sampling interval.

\textbf{Residual (shortcut) structure and interpretability.}
To explicitly incorporate a residual structure, the local state component in the input is passed through an identity shortcut
and added to obtain the one-step prediction.
Let $\hat{I}_i\in\mathbb{R}^{d\times(d(1+|Z_i|)+p+1)}$ be a linear selection matrix whose block form is defined by
\begin{equation}
\hat{I}_i = [I_d,\, 0_{d\times(d|Z_i|+p+1)}].
\label{eq:selector_matrix}
\end{equation}
This shortcut represents a persistence prior: over a short interval, the deviation state tends to change moderately,
and the network mainly learns the \emph{correction} caused by nonlinear rolling behavior and inter-stand coupling.

\textbf{Auxiliary branch for variable $\delta_n$.}
To improve robustness when the sampling interval length $\delta_n$ varies or becomes relatively large,
we introduce an auxiliary branch inside $\mathcal{N}_i$:
\begin{equation}
\mathcal{N}_i(X_{i,\text{in}};\Theta_i)\triangleq
\eta_i(X_{i,\text{in}};\Theta_{\eta_i}) + \mathcal{I}_i(X_{i,\text{in}};\theta_i),
\label{eq:aux_branch}
\end{equation}
where $\eta_i(\cdot)$ can be implemented by a lightweight feedforward network and
$\mathcal{I}_i(\cdot)$ denotes the remaining residual branch.
Conceptually, $\eta_i(\cdot)$ captures low-frequency/scale effects associated with interval length $\delta_n$,
while $\mathcal{I}_i(\cdot)$ learns the remaining coupling nonlinearities.
When $\eta_i(\cdot)\equiv 0$, the model reduces to the standard residual form.

Hence, the one-step prediction of the local deviation state at $t_{n+1}$ can be written as
\begin{equation}
X_{i,\text{out}} = \hat{I}_i X_{i,\text{in}} + \mathcal{N}_i(X_{i,\text{in}}; \Theta_i),
\label{333}
\end{equation}
where $X_{i,\text{out}}$ represents the predicted $\Delta x_i(t_{n+1})$.

\begin{remark}
The predictor in \eqref{333} admits a baseline-plus-correction form:
the shortcut term $\hat{I}_i X_{i,\mathrm{in}}$ propagates the current local deviation state $\Delta x_i(t_n)$,
while the residual network $\mathcal{N}_i(\cdot)$ learns the one-step correction.
This structure renders the model interpretable as a data-driven adjustment to a persistence prior,
with the correction capturing unmodeled nonlinearities and inter-stand coupling via $\Delta x_{Z_i}$
under varying operating conditions and varying sampling intervals.
\end{remark}

For the $j$-th one-step data sample, $j=1,\ldots,J$, we set
\begin{equation}
\begin{aligned}
X_{i,\text{in}}^{(j)} =
\big[
\Delta x_i^{(j)}(t_n),\ \Delta x_{Z_i}^{(j)}(t_n),\ 
\Gamma_{i,n}^{(j)},\ \delta_n^{(j)}
\big]^{\top}.
\end{aligned}
\label{eq:net_input_vector}
\end{equation}
The learning target is the one-step deviation-state change (residual)
\begin{equation}
\Delta r_i^{(j)}=\Delta x_i^{(j)}(t_{n+1})-\Delta x_i^{(j)}(t_n).
\label{eq:net_target_residual}
\end{equation}

%========================
\subsection{Training, Learned Model, and System Prediction (how prediction is used later)}
%========================

To suppress accumulation drift induced by long-horizon recursion and to improve long-term predictive stability,
we train the forward predictor jointly with an auxiliary backward residual model
and impose a multi-step reciprocal-consistency regularization over a $K$-step segment.

\textbf{Backward residual model.}
In addition to the forward residual predictor, we construct a backward residual network for subsystem $i$,
\begin{equation}
\mathcal{B}_i:\mathbb{R}^{d(1+|Z_i|)+p+1}\rightarrow\mathbb{R}^{d},
\end{equation}
parameterized by $\bar{\Theta}_i$. For the backward step associated with the interval $[t_n,t_{n+1}]$, we define
\begin{equation}
\begin{aligned}
X_{i,\mathrm{in}}^{b}
&=
\big[
\Delta x_i(t_{n+1}),\ \Delta x_{Z_i}(t_{n+1}),\
\Gamma_{i,n},\ \delta_n
\big]^{\top},\\
X_{i,\mathrm{out}}^{b}
&=
\hat{I}_i X_{i,\mathrm{in}}^{b} + \mathcal{B}_i(X_{i,\mathrm{in}}^{b};\bar{\Theta}_i),
\end{aligned}
\label{eq:backward_model_def}
\end{equation}
where $X_{i,\mathrm{out}}^{b}$ represents the backward estimate of $\Delta x_i(t_n)$.
Accordingly, the supervised backward residual target is
\begin{equation}
\Delta r_i^{b}=\Delta x_i(t_n)-\Delta x_i(t_{n+1}).
\end{equation}

\textbf{Forward rollout over a $K$-step segment.}
Given a segment sample $\mathcal{W}_{i,n}\in S_i^{(K)}$, we initialize the forward rollout by
\begin{equation}
\Delta \hat{x}_i(t_n)=\Delta x_i(t_n),
\end{equation}
and apply the forward predictor recursively for $K$ steps:
\begin{equation}
\begin{aligned}
\Delta \hat{x}_i(t_{n+s+1})
&=
\Delta \hat{x}_i(t_{n+s})
+
\mathcal{N}_i\!\Big(
\Delta \hat{x}_i(t_{n+s}),\,\Delta \hat{x}_{Z_i}(t_{n+s}),\,
\Gamma_{i,n+s},\,\delta_{n+s};\,\Theta_i
\Big), \\
&\qquad s=0,\ldots,K-1.
\end{aligned}
\label{eq:forward_rollout}
\end{equation}

\textbf{Backward rollout and reciprocal consistency.}
After obtaining the terminal forward prediction, we set the terminal condition for the backward rollout as
\begin{equation}
\Delta \bar{x}_i(t_{n+K})=\Delta \hat{x}_i(t_{n+K}),
\end{equation}
and roll back using $\mathcal{B}_i$ along the same segment:
\begin{equation}
\begin{aligned}
\Delta \bar{x}_i(t_{n+s})
&=
\hat{I}_i X_{i,\mathrm{in}}^{b}(t_{n+s})
+
\mathcal{B}_i\!\Big(X_{i,\mathrm{in}}^{b}(t_{n+s});\,\bar{\Theta}_i\Big),
\quad s=K-1,\ldots,0,
\end{aligned}
\label{eq:backward_rollout}
\end{equation}
where the backward input at time $t_{n+s}$ is
\begin{equation}
X_{i,\mathrm{in}}^{b}(t_{n+s})=
\big[
\Delta \bar{x}_i(t_{n+s+1}),\ \Delta \hat{x}_{Z_i}(t_{n+s+1}),\
\Gamma_{i,n+s},\ \delta_{n+s}
\big]^{\top},
\end{equation}
and $\Delta \hat{x}_{Z_i}(\cdot)$ is obtained from the same forward rollout.

With the forward and backward trajectories available on the same segment, we define the multi-step reciprocal prediction error
\begin{equation}
E_i(t_n)
=
\sum_{s=0}^{K}
\left\|
\Delta \hat{x}_i(t_{n+s})-\Delta \bar{x}_i(t_{n+s})
\right\|^2.
\end{equation}

\textbf{Training objectives (what each term enforces).}
We train the forward and backward networks jointly by minimizing the following objective terms:
\begin{equation}
\begin{aligned}
L_{\mathrm{1step}}(\Theta_i)
&= \frac{1}{J_K}\sum_{j=1}^{J_K}\frac{1}{K}\sum_{s=0}^{K-1}
\Big\|
\big(\Delta x_i^{(j)}(t_{n+s+1})-\Delta x_i^{(j)}(t_{n+s})\big)
-\mathcal{N}_i\!\left(
X_{i,\mathrm{in}}^{(j)}(t_{n+s});\Theta_i
\right)
\Big\|^2,\\[2mm]
L_{\mathrm{bwd}}(\bar{\Theta}_i)
&= \frac{1}{J_K}\sum_{j=1}^{J_K}\frac{1}{K}\sum_{s=0}^{K-1}
\Big\|
\big(\Delta x_i^{(j)}(t_{n+s})-\Delta x_i^{(j)}(t_{n+s+1})\big)
-\mathcal{B}_i\!\left(
X_{i,\mathrm{in}}^{b\,(j)}(t_{n+s});\bar{\Theta}_i
\right)
\Big\|^2,\\[2mm]
L_{\mathrm{msrp}}(\Theta_i,\bar{\Theta}_i)
&= \frac{1}{J_K}\sum_{j=1}^{J_K} E_i^{(j)}(t_n),\\[2mm]
L_{\mathrm{roll}}(\Theta_i)
&= \frac{1}{J_K}\sum_{j=1}^{J_K}\sum_{s=1}^{K}
\Big\|
\Delta x_i^{(j)}(t_{n+s})-\Delta \hat{x}_i^{(j)}(t_{n+s})
\Big\|^2.
\end{aligned}
\label{eq:training_losses}
\end{equation}
Here, $L_{\mathrm{1step}}$ enforces one-step accuracy, $L_{\mathrm{roll}}$ directly suppresses long-horizon error accumulation (drift),
and $L_{\mathrm{msrp}}$ regularizes the learned dynamics by reciprocal consistency between forward and backward rollouts.
In implementation, these terms can be combined by a weighted sum (weights chosen by validation):
\begin{equation}
L_{\mathrm{total}}=\lambda_1 L_{\mathrm{1step}}+\lambda_2 L_{\mathrm{roll}}+\lambda_3 L_{\mathrm{msrp}}+\lambda_4 L_{\mathrm{bwd}}.
\end{equation}

After training, we obtain the forward model
\begin{equation}
X_{i,\mathrm{out}} = \hat{I}_i X_{i,\mathrm{in}} + \mathcal{N}_i(X_{i,\mathrm{in}}; \Theta_i^*).
\end{equation}
For system prediction on a local interval $[t_n,t_{n+1}]$, we define the input vector as
\begin{equation}
X_{i,\mathrm{in}}
=
\big[
\Delta x_i(t_n),\ \Delta x_{Z_i}(t_n),\ \Gamma_{i,n},\ \delta_n
\big]^{\top},
\end{equation}
and perform the one-step prediction
\begin{equation}
\Delta \hat{x}_i(t_{n+1})
=
\Delta x_i(t_n)
+
\mathcal{N}_i\!\Big(
\Delta x_i(t_n),\,\Delta x_{Z_i}(t_n),\,
\Gamma_{i,n},\,\delta_n;\,\Theta_i^*
\Big).
\label{eq:one_step_predict}
\end{equation}
By applying this predictor recursively, we obtain a network model that predicts the system trajectory over long horizons.
Finally, the network parameters are optimized using the Adam optimizer:
\begin{equation}
\Theta_{i,t+1} = \Theta_{i,t} - \eta \frac{\hat{m}_{i,t}}{\sqrt{\hat{v}_{i,t}} + \varepsilon},
\end{equation}
where $\Theta_{i,t}$ denotes the current parameters, $\Theta_{i,t+1}$ the updated parameters, $\eta$ the learning rate,
$\hat{m}_{i,t}$ and $\hat{v}_{i,t}$ the bias-corrected first and second moment estimates, and $\varepsilon$ a small constant for numerical stability.
Figure~\ref{fig:rnn_logic} illustrates the overall structure.

\begin{figure}[htbp]
  \centering
  \includegraphics[scale=0.85]{picture/x6.pdf}
  \caption{Logic diagram of the residual neural network.}
  \label{fig:rnn_logic}
\end{figure}

%========================
\section{Nash Equilibrium-Based RNE-DMPC (five-stand, thickness--tension, clear objective and solution)}
%========================

The five-stand tandem cold rolling system is strongly coupled through inter-stand tension propagation.
As a result, changes in operating conditions or control actions (roll gap and speed) at one stand can affect both upstream and downstream stands,
making centralized online optimization over all decision variables computationally demanding.

To mitigate this issue, we decompose the global predictive-control problem into $N=5$ local subproblems associated with individual stands.
Each local controller optimizes its own decision variables while accounting for coupling via limited information exchange with neighboring controllers.
Motivated by game-theoretic coordination \citep{rawlings2008coordinating}, we formulate the distributed coordination process as a Nash-equilibrium-seeking iteration.

Based on the trained residual neural network surrogate model, we construct a Nash-equilibrium-based distributed model predictive control method (RNE-DMPC)
to achieve coordinated thickness--tension regulation/tracking. The overall control structure is shown in Figure~\ref{4}.

\begin{figure*}[htbp]
  \centering
  \includegraphics[width=\linewidth]{picture/x2.pdf}
  \caption{Schematic diagram of the control architecture for a tandem cold rolling mill.}\label{4}
\end{figure*}

\subsection{Control objective, references, and coupling conflict (what is optimized and why there is conflict)}
\label{subsec:objective_reference_conflict}

\textbf{Tracking variables.}
The controlled outputs are the exit thicknesses $\{h_i\}_{i=1}^5$ and inter-stand tensions $\{T_i\}_{i=1}^4$.
In deviation coordinates, the controller aims to drive $\Delta h_i(t)\to 0$ and $\Delta T_i(t)\to 0$ (regulation),
or track given time-varying references $h_i^{\mathrm{ref}}(t)$ and $T_i^{\mathrm{ref}}(t)$.

\textbf{Reference signal definition (clear meaning of $\Delta x_{\mathrm{ref}}$).}
For each prediction step $s$, define the local deviation reference vector of stand $i$ as
\begin{equation}
\Delta x_{i,\mathrm{ref}}(t_{n+s}) \triangleq
\begin{bmatrix}
h_i^{\mathrm{ref}}(t_{n+s})-h_i^{\mathrm{ref}}(t_{n+s})\\
T_{i-1}^{\mathrm{ref}}(t_{n+s})-T_{i-1}^{\mathrm{ref}}(t_{n+s})\\
T_{i}^{\mathrm{ref}}(t_{n+s})-T_{i}^{\mathrm{ref}}(t_{n+s})
\end{bmatrix}
=
\mathbf{0},
\end{equation}
when performing regulation around the references (i.e., we penalize deviation errors).
Equivalently, if one prefers to use absolute states, one may set the tracking error as
$\hat e_{x,i}(t_{n+s})=\hat x_i(t_{n+s})-x_i^{\mathrm{ref}}(t_{n+s})$.
In this paper, we keep the deviation-state formulation and directly penalize $\Delta\hat x_i(t_{n+s})$.
Thus, in the cost below, $\Delta x_{\mathrm{ref}}(t_{n+s})$ should be interpreted as the desired deviation (typically zero).
For boundary terms ($T_0,T_5$), their references are set to zero and their weights can be set to zero.

\textbf{Why conflict occurs (explicit coupling).}
Inter-stand tensions are shared coupling variables: the tension $T_i$ depends on the speeds of both stand $i$ and $i+1$ and the strip transport,
therefore it is affected by decisions from \emph{two} neighboring controllers.
At the same time, each controller also tries to achieve its \emph{local} thickness target via its roll gap.
Consequently, a speed/gap action that improves local thickness may deteriorate a shared tension, and vice versa.
This creates an intrinsic multi-agent conflict, motivating Nash-equilibrium coordination.

\subsection{Subsystem prediction and local optimization (decision variables, model constraint, cost, constraints, and solver)}
\label{subsec:local_mpc}

Each rolling stand is treated as a subsystem.
The coupled subsystem dynamics can be conceptually written as
\begin{equation}
\Delta x_i(t_{n+1})
=
f_i\!\big(\Delta x_i(t_n),u_i(t_n)\big)
+
\sum_{k \in Z_i} g_{ik}\!\big(\Delta x_k(t_n),u_k(t_n)\big),
\end{equation}
where $Z_i$ is given by \eqref{eq:neighbor_set} and $g_{ik}$ characterizes the coupling effect (mainly through tensions).
Instead of using an explicit first-principles model of $f_i,g_{ik}$, we use the trained neural-network surrogate to predict the deviation-state evolution.

\textbf{Local polynomial parameterization of control increments (two inputs).}
Over each local interval $[t_{n+s},t_{n+s+1}]$ with length $\delta_{n+s}$, the control increment trajectory of subsystem $i$
is parameterized by the vector quadratic polynomial \eqref{eq:du_poly_vector}:
\begin{equation}
\Delta u_{i,n+s}(\tau;\Gamma_{i,n+s})
=
\Gamma_{i,n+s,0}
+\Gamma_{i,n+s,1}\tau
+\Gamma_{i,n+s,2}\tau^2,
\qquad \tau \in [0,\delta_{n+s}],
\label{eq:du_poly_for_mpc}
\end{equation}
where $\Gamma_{i,n+s}\in\mathbb{R}^{p}$ with $p=6$ in \eqref{eq:Gamma_dim}.

\textbf{Neural-network-based prediction (model constraint for MPC).}
Using the trained residual neural network surrogate, subsystem $i$ predicts its one-step deviation state by
\begin{equation}
\begin{aligned}
\Delta \hat{x}_i(t_{n+s+1})
&=
\Delta \hat{x}_i(t_{n+s})
+
\mathcal{N}_i\!\Big(
\Delta \hat{x}_i(t_{n+s}),\,
\Delta \hat{x}_{Z_i}(t_{n+s}),\,
\Gamma_{i,n+s},\,
\delta_{n+s};\,
\Theta_i^*
\Big), \\
&\qquad s=0,\ldots,N_p-1,
\end{aligned}
\label{eq:nn_predict_in_mpc}
\end{equation}
where $\Delta \hat{x}_{Z_i}(t_{n+s})$ is obtained from the latest communicated neighbor predictions (Nash iteration described later).
Equation \eqref{eq:nn_predict_in_mpc} is the key link between \textbf{prediction} and \textbf{control}:
for any candidate decision variables $\{\Gamma_{i,n+s}\}$, the MPC optimizer can roll out \eqref{eq:nn_predict_in_mpc} to evaluate the predicted thickness/tension deviations.

\textbf{Decision variables.}
Optimize the local parameter sequence over the control horizon $N_c$:
\begin{equation}
\begin{aligned}
\mathbf{\Gamma}_i(t_n)
&=
\big[
\Gamma_{i,n}^\top,\,
\Gamma_{i,n+1}^\top,\,
\ldots,\,
\Gamma_{i,n+N_c-1}^\top
\big]^\top
\in \mathbb{R}^{pN_c}.
\end{aligned}
\end{equation}

\textbf{Local objective (explicit thickness--tension meaning).}
Let $\Delta \hat{x}_i(t_{n+s})=[\Delta \hat h_i(t_{n+s}),\,\Delta \widehat T_{i-1}(t_{n+s}),\,\Delta \widehat T_{i}(t_{n+s})]^\top$.
The local cost function of subsystem $i$ is defined as
\begin{equation}
\begin{aligned}
J_i
&=
\sum_{s=1}^{N_p}
\big\|
\Delta \hat{x}_i(t_{n+s}) - \Delta x_{i,\mathrm{ref}}(t_{n+s})
\big\|_{Q_i}^2
+
\sum_{s=0}^{N_c-1}
\big\|
\Gamma_{i,n+s}
\big\|_{R_i}^2.
\end{aligned}
\label{eq:local_cost}
\end{equation}
In the deviation-state regulation setting, $\Delta x_{i,\mathrm{ref}}(t_{n+s})\equiv 0$, so the first term penalizes predicted thickness and tension deviations.
The weighting matrix $Q_i$ is chosen to reflect the relative importance between thickness and tension regulation
(e.g., a larger weight on the first component for strict thickness control, while still penalizing adjacent tension deviations),
and $R_i$ penalizes the control-increment trajectory parameters to ensure smooth actuation.

\textbf{Constraints (how $\Gamma$ enforces bounds over the whole interval).}
Typical constraints include bounds on absolute inputs and increment trajectories:
\begin{align}
u_{i,\min} &\le u_i(t_{n+s}) \le u_{i,\max},
\label{eq:u_abs_bound}\\
\Delta u_{i,\min}
&\le
\Delta u_{i,n+s}(\tau;\Gamma_{i,n+s})
\le
\Delta u_{i,\max}, \notag \\
&\hspace{3.6cm}
\forall\tau\in[0,\delta_{n+s}].
\label{eq:du_traj_bound}
\end{align}
Here, $u_{i,\min},u_{i,\max}\in\mathbb{R}^{2}$ specify bounds on roll gap and speed (component-wise),
and $\Delta u_{i,\min},\Delta u_{i,\max}\in\mathbb{R}^{2}$ specify allowable increment bounds.
In practice, for the quadratic trajectory \eqref{eq:du_poly_for_mpc}, the interval-wise bound \eqref{eq:du_traj_bound} can be checked
by evaluating $\tau=0$, $\tau=\delta_{n+s}$, and the quadratic extremum $\tau^\star=-\Gamma_{i,n+s,1}\oslash(2\Gamma_{i,n+s,2})$
(component-wise) whenever $\tau^\star\in[0,\delta_{n+s}]$.

\textbf{Absolute-input update using interval average (consistency with execution).}
To enforce the absolute-input constraints consistently within the prediction horizon, we update the absolute input using the interval average:
\begin{equation}
\begin{aligned}
\Delta u_i(t_{n+s})
&=
\frac{1}{\delta_{n+s}}\int_0^{\delta_{n+s}}\Delta u_{i,n+s}(\tau;\Gamma_{i,n+s})\,d\tau \\
&=
\Gamma_{i,n+s,0}
+\Gamma_{i,n+s,1}\frac{\delta_{n+s}}{2}
+\Gamma_{i,n+s,2}\frac{\delta_{n+s}^2}{3},
\end{aligned}
\label{eq:du_avg_in_mpc}
\end{equation}
and propagate
\begin{equation}
\begin{aligned}
u_i(t_n) &= u_i(t_{n-1}) + \Delta u_i(t_n), \\
u_i(t_{n+s}) &= u_i(t_{n+s-1}) + \Delta u_i(t_{n+s}), \qquad s=1,\ldots,N_p-1.
\end{aligned}
\label{eq:u_update_in_mpc}
\end{equation}
This ensures smooth input evolution and avoids abrupt actuator changes.

\textbf{Local optimization problem (what is solved at each Nash iteration).}
At Nash-iteration index $l$, subsystem $i$ solves the following nonlinear optimization problem:
\begin{equation}
\mathbf{\Gamma}_i^{(l)}=\arg\min_{\mathbf{\Gamma}_i}\; J_i
\quad \text{s.t.}\quad
\eqref{eq:nn_predict_in_mpc},\ \eqref{eq:u_abs_bound}\text{--}\eqref{eq:u_update_in_mpc}.
\label{eq:local_nlp}
\end{equation}
Because the neural surrogate $\mathcal N_i(\cdot)$ is differentiable, \eqref{eq:local_nlp} is a differentiable nonlinear program (NLP).
It can be solved by standard gradient-based NLP solvers (e.g., SQP or interior-point methods) using automatic differentiation to compute gradients.

\subsection{Nash equilibrium coordination iteration (how conflict is resolved and how the final solution is obtained)}
\label{subsec:nash_coordination}

The Nash equilibrium is computed via distributed best-response iterations. Each controller repeatedly computes its best response to the latest neighbor strategies
and exchanges predicted trajectories and decision variables through a communication module.
This resolves coupling conflicts by seeking a strategy profile in which no single subsystem can unilaterally reduce its own objective \eqref{eq:local_cost}
given the other subsystems' strategies.

The distributed best-response iteration is summarized in Table~\ref{tab:nash_iter_en}.

\begin{table}[t]
\centering
\small
\renewcommand{\arraystretch}{1.12}
\setlength{\tabcolsep}{3.5pt}
\caption{Distributed Nash best-response iteration for RNE-DMPC (five-stand).}
\label{tab:nash_iter_en}
\begin{tabularx}{\linewidth}{>{\centering\arraybackslash}p{0.11\linewidth} X}
\toprule
\textbf{Step} & \textbf{Description} \\
\midrule
A &
Initialize $l=1$ and initialize $\mathbf{\Gamma}_i^{(0)}$ for all stands (e.g., warm-start from previous time step). \\

B &
Using the surrogate predictor \eqref{eq:nn_predict_in_mpc}, compute $\Delta \hat{x}_i^{(l)}(t_{n+s})$ for $s=1,\ldots,N_p$ \\
&
given $\mathbf{\Gamma}_i^{(l-1)}$ and the latest neighbor predictions
$\Delta \hat{x}_{Z_i}^{(l-1)}(t_{n+s})$. \\

C & Solve the local NLP \eqref{eq:local_nlp} to update $\mathbf{\Gamma}_i^{(l)}$ (best response). \\

D &
Broadcast $\mathbf{\Gamma}_i^{(l)}$ and predicted trajectories
$\Delta \hat{x}_i^{(l)}(t_{n+s})$ to the communication system. \\

E &
Update neighbor predictions $\Delta \hat{x}_{Z_i}^{(l)}(t_{n+s})$ using received information; re-generate predictions if needed. \\

F & Compute the maximum relative change $\varsigma^{(l)}$ as the convergence metric. \\

G &
If $\varsigma^{(l)}\le \varsigma_{\mathrm{tol}}$, stop and set $\mathbf{\Gamma}_i^*=\mathbf{\Gamma}_i^{(l)}$; \\
& otherwise set $l\leftarrow l+1$ and repeat Steps B--F. \\
\bottomrule
\end{tabularx}
\end{table}

The convergence metric in Step~F is defined as
\begin{equation}
\begin{aligned}
\varsigma^{(l)}
&=
\max_i
\frac{\left\|
\mathbf{\Gamma}_i^{(l)}-\mathbf{\Gamma}_i^{(l-1)}
\right\|_2}{
\left\|
\mathbf{\Gamma}_i^{(l-1)}
\right\|_2+\epsilon},
\end{aligned}
\end{equation}
where $\epsilon>0$ avoids division by zero.

\textbf{Receding-horizon implementation (how control is applied).}
After the Nash iteration converges, only the first-interval parameters $\Gamma_{i,n}^*$ are applied.
The control increment trajectory over $[t_n,t_{n+1}]$ is
\begin{equation}
\Delta u_{i,n}(\tau)=\Delta u_{i,n}(\tau;\Gamma_{i,n}^*),
\quad \tau\in[0,\delta_n].
\end{equation}
To obtain the equivalent discrete-time control increment used for updating the absolute input, we use the interval average:
\begin{equation}
\begin{aligned}
\Delta u_i(t_n)
&=
\frac{1}{\delta_n}\int_0^{\delta_n}\Delta u_{i,n}(\tau)\,d\tau \\
&=
\Gamma_{i,n,0}^*
+
\Gamma_{i,n,1}^*\frac{\delta_n}{2}
+
\Gamma_{i,n,2}^*\frac{\delta_n^2}{3}.
\end{aligned}
\end{equation}
Then the absolute control input is updated by
\begin{equation}
u_i(t_n)=u_i(t_{n-1})+\Delta u_i(t_n),
\end{equation}
which ensures smooth input evolution and avoids abrupt actuator changes.

\textbf{How prediction and control are connected in closed loop.}
At each sampling instant $t_n$:
(i) measure/estimate the current deviation states $\Delta x_i(t_n)$ (thickness and tensions) and obtain neighbor information;
(ii) run the Nash best-response iteration using the neural predictor \eqref{eq:nn_predict_in_mpc} as the model constraint to solve \eqref{eq:local_nlp};
(iii) apply the first control-increment trajectory $\Delta u_{i,n}(\tau;\Gamma_{i,n}^*)$ and update $u_i(t_n)$ by interval averaging;
(iv) roll the horizon forward to $t_{n+1}$ and repeat.
In this way, the learned predictor directly provides the multi-step forecasts needed by MPC, and the optimization computes actuator trajectories
(roll gap and speed increments) that minimize thickness--tension deviations while coordinating coupling conflicts via Nash equilibrium.

The overall control flow chart is shown in Fig.~\ref{liu}.

\begin{figure}[htbp]
  \centering
  \includegraphics[width=\linewidth]{picture/x5.pdf}
  \caption{The overall system control flow chart}\label{liu}
\end{figure}