作者: KSkun

NTP 网络时间协议

NTP 网络时间协议

By KSkun, 2020/7/20 NTP(Network Time P 

QUIC 快速 UDP 互联网连接协议

QUIC 快速 UDP 互联网连接协议

KSkun 2020/6/12

2013 年,Google 提出了快速 UDP 互联网连接协议(QUIC),实现了一种可以替代 TLS/TCP 模式的传输协议,并且将其应用在 Google 搜索、YouTuBe、Chrome 浏览器等产品上。在这些产品上收集到的数据表明,QUIC 相较于传统的 TLS/TCP 模式具有较大的性能优势。2015 年,QUIC 被提交至 IETF,随后在 2018 年被确定为 HTTP/3 标准中的内容,正式纳入 HTTP 家族。

为什么 Google 要开发 QUIC?为什么选择了 UDP?QUIC 如何优化性能?QUIC 效果如何?接下来,我们将一一解答这些问题,并阐释 QUIC 的基本内容。

现在的 HTTPS 如何工作?

要了解 QUIC 的开发背景,首先需要知道今天的 HTTP 以什么形式工作。因此,我们将介绍 HTTPS 协议栈及其存在的问题,从而引出 QUIC 的优化。

传统 HTTPS 协议栈

传统的 HTTPS 协议栈包含以下三部分:

  • 应用层:HTTP/2,提供 HTTP 应用功能,由客户端(浏览器)在用户空间实现;
  • 安全层:TLS,提供安全的传输链路,由客户端(浏览器)在用户空间实现;
  • 传输层:TCP,提供可靠传输方式,由操作系统在内核中实现。

在最近的几年中,上述几种协议均进行了许多改进,例如:

  • HTTP/2 解决了 HTTP/1.1 的性能问题,对 HTTP 头进行压缩,引入了 Server Push(服务器主动推送流,而不是等客户端发送请求后再响应,减少了网页加载多次请求资源的性能问题),在单个 TCP 连接里实现多路复用请求等;
  • TLS 1.3 也优化了 TLS 1.2 存在的一些问题,包括废弃了过时的加密算法,支持 0-RTT 传输(优化握手时间)等;
  • TCP 中,CUBIC、BBR 等优秀的拥塞控制算法提高了网络资源的利用率,优化了传输效率。

容易发现,HTTP 与 TLS 都进行了较大规模的改动,而 TCP 并没有对其协议本身进行修订,这一现象是存在原因的。

「僵化」与「耦合」

在 QUIC 的论文中,Google 使用了僵化(ossification)一词来形容现在的互联网中间设备。这是因为,许多设备供应商在设备中预置了一些规则,如防火墙会禁止已有规则之外的传输,NAT 会以一定规则重写传输头等。如果改动 TCP 协议,这些设备必须更新规则才能支持变动,设备更新的不及时则会影响修订的推进。

此外,由于 TCP 通常实现在操作系统内核中,这种 TCP 与操作系统的耦合(coupling)导致如果修改 TCP 协议,则必须升级操作系统内核,而内核的改动与升级的下发的时间周期很长,也影响了修订的推进。

协议本身存在的问题

由于 TLS 是在 TCP 建立连接后再进行握手,除了 TCP 建立连接时的 1-RTT 时间开销,TLS 还增加了一个 2-RTT 的握手环节。Web 连接大多都是短暂的连接,频繁进行不必要的重复握手对短传输性能存在影响。

QUIC8 - QUIC 快速 UDP 互联网连接协议
图 1:TLS 握手的流程

在 HTTP/1.1 和 HTTP/2 中,协议限制了最大 TCP 连接数量,而 HTTP/2 使用了多路复用同一 TCP 连接的技术,这带来了新的问题:行头阻塞(head-of-line blocking)。当复用同一 TCP 连接时,由于要求响应按顺序传输,第一个响应丢包重传时,会阻塞后面耗时较短的响应及时传输。

QUIC9 - QUIC 快速 UDP 互联网连接协议
图 2:行头阻塞的示意图

QUIC 是什么?

广泛使用的 TLS/TCP 仍存在很大的改进空间,但由于种种困难无法推行。为了彻底解决 TCP 带来的麻烦,Google 选择了 UDP,继而创造出了本文的主角——QUIC。

QUIC(Quick UDP Internet Connection,快速 UDP 互联网连接协议)是一种以 UDP 为底层传输协议,支持加密、多路复用,工作在用户空间的的低延迟传输协议。

QUIC1 1024x512 - QUIC 快速 UDP 互联网连接协议
图 3:QUIC 相对于传统 HTTPS 协议栈的地位

如图所示,QUIC 取代了传统 HTTPS 协议栈中部分 TCP、TLS 与部分 HTTP/2 的位置,向下使用 UDP 进行传输,向上提供 HTTP 接口。QUIC 具有与 HTTPS 相同的接口与安全特性,且低延迟、高效率、能够进行快速的迭代更新,有效解决了上面所提到了诸多问题。

接下来,我们通过 Google 提供的实验数据了解 QUIC 优化传输效率的效果。Google 通过握手延迟、搜索延迟(在 Google 搜索引擎上实验)、视频延迟、重缓冲率(视频暂停缓冲的时间占总时间的比例,皆在 YouTuBe 上实验)来考察 QUIC 对优化不同应用场景性能的效果。

QUIC2 1024x384 - QUIC 快速 UDP 互联网连接协议
图 4:在不同 RTT 链路情况中,各种协议的握手延迟
QUIC3 1024x384 - QUIC 快速 UDP 互联网连接协议
图 5:在不同 RTT 链路情况中,各种协议的搜索延迟
QUIC4 1024x384 - QUIC 快速 UDP 互联网连接协议
图 6:在不同 RTT 链路情况中,各种协议的视频延迟
QUIC5 1024x384 - QUIC 快速 UDP 互联网连接协议
图 7:在不同 RTT 链路情况中,各种协议的视频重缓冲率

可以发现,在各种情况下,QUIC 的性能不同程度优于 TCP/TLS 的表现。

QUIC 如何工作?

QUIC 是如何实现优秀的性能的?要解答这个问题,就必须理解 QUIC 的工作机制,并将其与 TLS/TCP 模式对比,从而体现 QUIC 的优势。接下来,我们将简要介绍 QUIC 的几个核心机制。

建立连接

QUIC6 - QUIC 快速 UDP 互联网连接协议
图 8:三种不同的 QUIC 的握手情形

上图描述了三种 QUIC 的握手情形,接下来分情况研究这三种情形。

初始握手:第一次握手之前,客户端并不知道任何有关服务器的信息,因此客户端发送一个不完整的客户端问候(CHLO)信息,服务器回应一个拒绝(REJ)信息。服务器回应的 REJ 中包含了服务器的配置细节(包含服务器长期公钥)、证书链、签名、源地址令牌。源地址令牌包含了客户端的公共 IP 地址与时间戳,客户端在接下来的信息中包含这一令牌来代表自己的身份。客户端通过签名来验证服务器发送的 REJ 信息,并将使用这些信息发送完整的 CHLO(包含客户端的临时公钥)。

最终握手(或重复握手):初始握手使客户端获得了与服务器建立连接的所有细节,通过发送完整的 CHLO 信息与服务器成功建立连接。此时客户端使用的密钥为初始密钥,包含客户端的临时私钥和服务器的长期公钥。同时,为了达到 0-RTT 的握手延迟,客户端会同时开始使用已知信息发送数据,而非等待服务器回应后再发送。

如果握手成功,服务器回复一个服务器问候(SHLO)信息,这个信息使用初始密钥加密后发送,包含了服务器的临时公钥。客户端收到 SHLO 后,将切换使用服务器的临时公钥加密。这样一来,除了客户端最初的数据使用初始密钥加密,其他数据都使用临时密钥加密,实现了数据的前向安全性(即使服务器长期密钥泄露,也无法使之前的数据加密失效)。

客户端会缓存服务器配置与源地址令牌,在下一次连接同一服务器时,直接以缓存的配置与服务器进行握手,实现 0-RTT 握手。

客户端缓存的配置或令牌有可能过期,此时服务器会直接发送 REJ 来拒绝握手,并附带最新的信息。客户端从 REJ 中提取最新信息,更新缓存后即可再发起握手。

版本协商:客户端会首先提议一个版本,并以该版本传输数据。如果服务器不支持客户端提议的版本,则发回一个强制版本协商数据,包含了服务器支持的所有版本。如果服务器的版本总是高于客户端版本,则版本协商是 0-RTT 的,因此鼓励服务器及时更新 QUIC 版本。为了防止降级攻击,生成密钥的时候将把客户端与服务器发送的版本信息也纳入参数之中,避免信息被修改。

多路复用

为了避免产生行头阻塞的问题,QUIC 支持在同一连接中同时传输多个流,丢包只影响到涉及到的流,而其他流仍然能正常传输而不被阻塞。

QUIC 流是可靠的双向字节流,单个流最高可以传输 264 字节数据。流非常轻量,因此发送少量数据时也可以使用新流。流通过流 ID 进行标识,为了避免 ID 冲突,客户端发起的流 ID 为奇数,服务器为偶数。流可以在发送数据时隐式发起,通过在最后一个流帧中设置 FIN 标记关闭。如果流不再需要,可以取消流而不断开 QUIC 连接。

QUIC 数据包的格式如下,一个数据包包括公共头与多个流帧,每个流帧具有自己的头与数据,数据包可以承载来自不同流的流帧。

QUIC7 1024x673 - QUIC 快速 UDP 互联网连接协议
图 9:QUIC 数据包的结构

QUIC 发送数据的速率是受限的,因此必须决定如何在多个流之间分配速率。在论文实现中,QUIC 使用 HTTP/2 流优先级来分配速率。

身份验证和加密

QUIC 数据包中未加密的部分(包头等)对路由或解密数据都是必要的。Flags 编码了连接 ID 与包编号的长度;连接 ID 可以被负载均衡使用来分配服务器或被服务器用于识别连接状态;版本号和多样化 nonce 只出现在较早的包中,服务器生成 nonce 为客户端生成密钥增加熵值。两方都以包编号作为每个包的 nonce,用于进行加密与解密。

在未加密的握手包(如版本协商信息)中包含的信息都被用于生成最终密钥,篡改这些包会导致生成的密钥不一致,从而导致连接失败。当服务器中不存在连接状态时,会发出重置包,通常路由更改或服务器重启会造成这种情况。由于服务器不知道连接的密钥,重置包只能不加密且不经验证。

丢包恢复

由于 TCP 重传时,ACK 包的编号与原始编号相同,ACK 接受者无法判断确认的是重传还是原始包,必须等待 1-RTT 才能确认这一信息。每个 QUIC 包都有不同的编号,重传包也是一样,因此不需要区分是原始还是重传包。在 QUIC 包中,包编号表示时间顺序,而流帧中的偏移量表示数据的顺序,这些信息可以更准确地检测丢包。

QUIC ACK 包含了接受数据与发送 ACK 之间的延迟,包含这一信息有助于更准确地估计 RTT,为 BBR 等拥塞控制算法提供帮助。

流量控制

QUIC 同时限制了单个流可以消耗的缓冲区大小与一个连接的缓冲区大小,避免了流与连接之间的行头阻塞。

QUIC 接收方会首先发送每个流中的缓冲窗口大小,并定期在流中发送窗口更新帧来扩大窗口大小限制。连接级窗口也用相似的方法控制。此外,实现还采用了类似 TCP 中的调整。

拥塞控制

QUIC 不依赖于某一特定拥塞控制算法,而是提供控制接口,支持多种算法。

连接迁移

QUIC 连接由一个 64 位的连接 ID 标识,该 ID 可以保存连接的状态。当客户端连接发生迁移时,可以通过连接 ID 来迁移连接的状态,避免 NAT 更改丢失状态的问题。

QUIC 发现

客户端并不知道服务器是否支持 QUIC,因此客户端第一次发送 HTTP 请求时先使用 TLS/TCP,而服务器在响应中附带一个 Alt-Svc 头来说明支持 QUIC。

在后序请求中,客户端争用 QUIC 和 TLS/TCP,但尝试 QUIC 略微超前 TLS/TCP。如果 QUIC 无法连接,则客户端之后使用 TLS/TCP。

QUIC 的局限性

尽管在实验中表现优秀,QUIC 仍然具有其局限性。Google 在其论文中提到了以下几点问题:

超前连接:应用程序通过提前进行握手来避免握手造成的延迟,这种情况下 QUIC 的 0-RTT 握手并不具有太大的优势。

高带宽、低延迟、低丢包率的网络:在这样的网络中使用 QUIC 几乎没有优势,反而有可能带来性能损失。

移动设备:QUIC 在移动设备中优势也不大,通常是因为移动设备应用针对其网络特点进行了优化。而 QUIC 需要较多的 CPU 性能,在移动设备上表现并不优。

UDP 限制:部分网络会限制或禁用 UDP 流量,QUIC 无法建立连接。

网络中间设备的僵化:QUIC 是一种快速迭代的协议,中间设备有时会附带对 QUIC 连接的限制规则,当 QUIC 更新时设备规则并未及时更新,从而影响到 QUIC 的正常连接。

CPU 占用较高:QUIC 的 CPU 使用率约是 TLS/TCP 的 3.5 倍。Google 通过优化加密算法等方式降低 CPU 使用率至 TLS/TCP 的 2 倍,但仍然较高。

结语

Google 对互联网技术发展做出的贡献巨大,相继创造了 BBR 拥塞控制算法、QUIC 协议等更高效的互联网技术,并推动它们成为现代互联网基础设施的一部分。

截至目前,QUIC 还未得到大规模的部署与使用。我们期待着 QUIC 能够在未来的互联网中充分发挥其性能优势,在互联网中创造更多可能性;同时也期待着更多优秀的互联网技术被创造、被接纳、被广泛使用,让世界更加触手可及。

参考资料

本文引用的图片部分来源于论文或互联网,引用来源均列入参考资料列表中,感谢原作者的分享。

BBR(瓶颈带宽和往返传播时间)拥塞控制算法

BBR(瓶颈带宽和往返传播时间)拥塞控制算法

为什么需要拥塞控制?

在网络上传输数据就像在水管中传输水流。假如我们要通过一根参数未知的水管输水,输得多了水管会爆裂,输得少了效率会降低,如何传输就是一个关键的问题。

  • 最开始的时候,我们不知道这根水管能承受多大的速率,我们可以监控水管的状况,并且逐渐开大阀门,直到发现状况异常;
  • 如果输水过程中,水管被砸了一下而变窄了(或者突然改道变宽了),以前的速率可能过大(或过小),为了适应新的参数,我们可能需要调整速率,比如把阀门关上重新逐渐开大。

也许上面的策略不是最好的策略,我们还可以制定更好的策略让调整的时间变短,从而使输水的效率增加。这种调整的策略就是拥塞控制算法

v2 6a3a7b0f85bf9a664e8f9a000d80a8e7 r - BBR(瓶颈带宽和往返传播时间)拥塞控制算法

根据上述场景,我们可以确定理想的拥塞控制算法满足的要求:

  • 能在尽量短时间内获得当前链路的传输能力;
  • 能及时适应链路的波动,确定较合适的传输速率。

以前我们是怎么进行拥塞控制的?

Tahoe 算法

  • 拥塞窗口:最大允许没有收到 ACK 确认的包的数量
  • 慢启动:每收到一个 ACK,拥塞窗口指数增长,直到丢包或达到慢启动阈值
  • 发生丢包时,重传丢包并重新慢启动,达到慢启动阈值后线性增长

Reno 算法

  • 发生丢包时,拥塞窗口减半并直接线性增长
v2 86c12ed3ddb2391c057f7d3598372880 720w - BBR(瓶颈带宽和往返传播时间)拥塞控制算法

New Reno 算法

  • 仅对重传部分的效率有优化

CUBIC 算法

  • 窗口大小是上次拥塞事件以来的时间的三次函数
  • 利用三次函数的凹部分和凸部分实现快速启动、平缓增长
  • 不依赖 ACK 包,与延迟无关
  • 适应于长胖网络:高带宽、高延迟

局限性

容易发现,上述几种算法都基于对丢包的判断和处理来进行拥塞控制。但是现在,丢包和拥塞的关系不那么紧密了,这些算法对于高延迟、高带宽、有一定丢包率的网络并不能达到最大传输效率。

怎么解决基于丢包的拥塞控制算法的问题?

别老盯着丢包!

既然丢包不是一种合适的判断方法,那么我们干脆不以丢包作为拥塞控制的依据,而是另寻其他信息。

Google 的研究人员对世界各地的大量 TCP 头部进行研究,确定了一种新的判断依据:瓶颈带宽(Bottleneck Bandwidth)和往返传播时间(Round-trip Propagation Time)。

  • 瓶颈链路:在网络中,一个连接的一段最慢的链路
  • 瓶颈带宽:瓶颈链路的带宽
  • 往返传播时间

考虑正常行驶的高速公路上突然发生了交通事故,这会导致局部通行速率降低,使得该处有车辆形成队列,从而使整条道路的通信速率降低。瓶颈带宽就是事故处的通行速率,而往返传播时间可以代表道路的长度。

瓶颈带宽、往返传播时间和流行为的关系

vanjacobson1 - BBR(瓶颈带宽和往返传播时间)拥塞控制算法

丢包往往发生在图中缓存限制的区域,因此传统的 TCP 拥塞控制算法事实上只能将传输情况控制在该区域的边界附近。在此处,网络的延迟较大,缓存队列较慢,并不是最优方案。

通过 BBR 算法,我们可以以瓶颈带宽、往返传播时间来控制传输,使得网络中不存在缓存队列。这相当于控制向水管注水的速率,使水管中最细的一段恰好填满且没有水堆积在近端。

计算出瓶颈带宽、往返传播时间

往返传播时间是链路的物理特性,可以视作不会改变,因此可以统计一段时间内的延迟最小值作为往返传播时间。

我们无法直接对瓶颈带宽进行计算,但可以用间接的方法。可以统计一个包从发出到收到 ACK 的时间间隔,在统计这段时间内的未确认包数,就可以计算出瓶颈带宽。

适应网络的参数变化

由于网络可能不断产生变化,瓶颈带宽、往返传播时间也应该进行实时的估计与调整。但当网络带宽增大时,上述过程并不会主动探测到带宽的增大。

BBR 使用 pacing 作为试探瓶颈带宽的方式,每次将窗口设置为比当前瓶颈带宽稍大的值,重新统计速率,即可判断当前状况处于临界点的哪一边。如果处于临界点右侧,则之后使用比瓶颈带宽稍小的值来消除网络中的队列。

vanjacobson2 - BBR(瓶颈带宽和往返传播时间)拥塞控制算法
vanjacobson3 - BBR(瓶颈带宽和往返传播时间)拥塞控制算法

连接的启动

BBR 在连接启动时,使用指数上升的规律增大发送速率,直到探测出瓶颈带宽。一旦找到瓶颈带宽,BBR 使用增大倍数的倒数消耗队列,之后保持上述稳定状态。

vanjacobson4 - BBR(瓶颈带宽和往返传播时间)拥塞控制算法

BBR 可以主动消耗网络中的队列,而传统的拥塞控制算法无法消耗网络中的队列,这是 BBR 的优势。

vanjacobson5 - BBR(瓶颈带宽和往返传播时间)拥塞控制算法

公平分配瓶颈链路

在新的连接建立之前,可能已有一些连接占满了瓶颈链路,BBR 需要平均分配链路,使得新建立的连接也可以分配到传输资源。

当一段时间内往返传播时间没有改变时,BBR 主动将窗口减小,使得网络中队列被消耗,其他连接得到了更小的往返传播时间,使得其他连接也主动减小窗口,之后同时开始回到原状态,从而实现资源的重新公平分配。

vanjacobson6 - BBR(瓶颈带宽和往返传播时间)拥塞控制算法

参考资料

Goodbye, 2019

Goodbye, 2019

这是一篇 KSkun &#3034 

[CSP-S2 2019]Emiya 家今天的饭 题解

[CSP-S2 2019]Emiya 家今天的饭 题解

题目地址:洛谷:P5664 Emiya 家今天的饭 – 洛谷 | 计算机科学教育新生态

题目描述

Emiya 是个擅长做菜的高中生,他共掌握 $n$ 种烹饪方法,且会使用 $m$ 种主要食材做菜。为了方便叙述,我们对烹饪方法从 $1 \sim n$ 编号,对主要食材从 $1 \sim m$ 编号。

Emiya 做的每道菜都将使用恰好一种烹饪方法与恰好一种主要食材。更具体地,Emiya 会做 $a_{i,j}$​ 道不同的使用烹饪方法 $i$ 和主要食材 $j$ 的菜($1 \leq i \leq n, 1 \leq j \leq m$),这也意味着 Emiya 总共会做 $\sum\limits_{i=1}^{n} \sum\limits_{j=1}^{m} a_{i,j}$​ 道不同的菜。

Emiya 今天要准备一桌饭招待 Yazid 和 Rin 这对好朋友,然而三个人对菜的搭配有不同的要求,更具体地,对于一种包含 $k$ 道菜的搭配方案而言:

  • Emiya 不会让大家饿肚子,所以将做至少一道菜,即 $k \geq 1$​
  • Rin 希望品尝不同烹饪方法做出的菜,因此她要求每道菜的烹饪方法互不相同
  • Yazid 不希望品尝太多同一食材做出的菜,因此他要求每种主要食材至多在一半的菜(即 $​\lfloor \frac{k}{2} \rfloor$​ 道菜)中被使用
    • 这里的 $​\lfloor x \rfloor$​ 为下取整函数,表示不超过 $​x$​ 的最大整数。

这些要求难不倒 Emiya,但他想知道共有多少种不同的符合要求的搭配方案。两种方案不同,当且仅当存在至少一道菜在一种方案中出现,而不在另一种方案中出现。

Emiya 找到了你,请你帮他计算,你只需要告诉他符合所有要求的搭配方案数对质数 $​998,244,353$​ 取模的结果。

题意简述

一共有 $n$ 种烹饪方法和 $m$ 种主要食材,一道菜仅对应一种烹饪方法和一种主要食材。使用第 $i$ 种烹饪方法和第 $j$ 种主要食材,Emiya 可以做出 $a_{i,j}$ 道不同的菜。

求满足以下条件的菜的集合数:

  • 非空;
  • 每道菜的烹饪方法互不相同;
  • 集合中每种主要食材的菜数不超过集合大小的一半(向下取整)。

输入输出格式

输入格式:

第 1 行两个用单个空格隔开的整数 $n,m$。

第 2 行至第 $n + 1$ 行,每行 $m$ 个用单个空格隔开的整数,其中第 $i + 1$ 行的 $m$ 个数依次为 $a_{i,1}, a_{i,2}, \cdots, a_{i,m}$​。

输出格式:

仅一行一个整数,表示所求方案数对 $998,244,353$ 取模的结果。

输入输出样例

样例 #1

输入样例#1:

2 3
1 0 1
0 1 1

输出样例#1:

3

样例 #2

输入样例#2:

3 3
1 2 3
4 5 0
6 0 0

输出样例#2:

190

样例 #3

输入样例#3:

5 5
1 0 0 1 1
0 1 0 1 0
1 1 1 1 0
1 0 1 0 1
0 1 1 0 1

输出样例#3:

742

数据范围

测试点编号$n=$$m=$$a_{i,j}<$
$1$$2$$2$$2$
$2$$2$$3$$2$
$3$$5$$2$$2$
$4$$5$$3$$2$
$5$$10$$2$$2$
$6$$10$$3$$2$
$7$$10$$2$$2$
$8$$10$$3$$1000$
$9 \sim 12$$40$$2$$1000$
$13 \sim 16$$40$$3$$1000$
$17 \sim 21$$40$$500$$1000$
$22 \sim 25$$100$$2000$$998,244,353$

对于所有测试点,保证 $1 \leq n \leq 100, 1 \leq m \leq 2000, 0 \leq a_{i,j} \lt 998,244,353$。

题解

参考资料:题解 P5664 【Emiya 家今天的饭【民间数据】】 – Caro23333 的博客 – 洛谷博客

居然在 NOIP 级别的比赛中见到了计数 DP,很难想象以后的 NOIP 对于数学知识的考察会变成什么样。那个只有一道很水的数论/组合签到题的时代已经结束了。

首先本题是组合数学/计数类型的问题,这类题通常有两种解决方式:推公式或 DP,本题中可以采用 DP 求解。

接下来考虑如何设计 DP 状态。如果正面解决,就必须对每种主要食材的出现次数做小于集合大小的一半的限制,这需要记录下每种主要食材的出现次数,不是很可行。因此可以采用正难则反的策略,容易发现不合法的情况中只会有一种主要食材超过集合大小的一半,处理这个限制仅需记下一种主要食材的出现次数,问题变得简单了。

所有方案数

首先从一个较好解决的问题入手,由于答案=所有方案数-不合法方案数,我们需要知道所有的方案有多少种。

定义 $g(i, j)$ 表示前 $i$ 种烹饪方式中,每种最多选出一道菜的情况下,当前集合中有 $j$ 道菜的方案数,则决策为当前处理的第 $i$ 种烹饪方式是否选出一道菜,记 $s(i)$ 代表第 $i$ 种烹饪方式一共可做出的菜数,即 $s(i) = \sum\limits_{j=1}^m a_{i, j}$ 转移方程如下:

$$ g(i,j) = \underset{不选}{g(i-1, j)} + \underset{从 s(i) 道菜中选 1 种}{g(i-1,j-1) \cdot s(i)} $$

初始值为 $g(0, 0)=1$(什么都没做,有 $1$ 种方案),所有的方案数为 $\sum\limits_{i=1}^n g(n, i)$,由于集合不能为空,$g(n, 0)$ 的值不能加入方案数中。

这一部分的复杂度为 $O(n^2)$。

不合法方案数

由于不合法的方案中,只有一种主要食材的数量会超过集合大小的一半,可以枚举超限的主要食材,而只记录下其他食材的数量之和,不分开记录。

定义 $f(in, j, k)$ 表示考虑前 $in$ 种烹饪方式,当前要超限的主要食材的出现次数为 $j$,其他主要食材的出现次数之和为 $k$。决策有三种:当前考虑的第 $in$ 种烹饪方式不选、选要超限的那种主要食材的菜或选其他主要食材的菜,记 $s(i)$ 代表第 $i$ 种烹饪方式一共可做出的菜数,且当前要超限的主要食材编号为 $im$,有转移方程:

$$ f(in, j, k) = \underset{不选}{f(in-1, j, k)} + \underset{选要超限的主要食材}{f(in-1, j-1,k) \cdot a_{in, im}} + \underset{选其他的主要食材}{f(in-1,j,k-1) \cdot [s(in)-a_{in,im}]} $$

将上述 DP 在以每一种主要食材为超限食材的情况下都计算一遍,对于 $im$ 这种食材为超限食材的情况,不合法的方案数为 $\sum\limits_{k < j} f(n, j, k)$。将这些不合法的方案数加起来就获得了总的不合法方案数了。但这种方式的复杂度是 $O(n^3m)$ 的,无法通过本题。

考虑优化上述方法,将不合法的方案满足的条件 $k < j$ 移项,得到表达式 $j-k > 0$,其实判断一个方案是否合法,仅需要判断 $j-k$ 的符号即可,我们考虑将之前状态中的 $j$ 和 $k$ 两个维度变为一维 $j-k$ 的值,这将会减少一层循环,使复杂度降至 $O(n^2m)$。

接下来考虑上面的改变会对转移方程造成什么样的影响,事实上,仅需将方程中的两维改为一维即可:

$$ f(in, j) = \underset{不选}{f(in-1, j)} + \underset{选要超限的主要食材}{f(in-1, j-1) \cdot a_{in, im}} + \underset{选其他的主要食材}{f(in-1,j+1) \cdot [s(in)-a_{in,im}]} $$

由于 $j-k$ 的值可以为负数,而数组下标不能为负数,因此下标应当处理为 $j-k+n$,以让下标的值非负。DP 的初始值为 $f(0, 0)=1$(什么都没做,其中第二维下标应处理为 $n$),对于每一种超限食材的情况,不合法的方案数为$\sum\limits_{j=1}^n f(n, j)$(第二维下标应处理为 $j + n$)。

最后,将所有方案数减去不合法方案数就得到了最终答案,总复杂度为 $O(n^2m)$。

代码

// Code by KSkun, 2019/12
#include <cstdio>
#include <cctype>
#include <cstring>

#include <algorithm>
#include <utility>

typedef long long LL;
typedef std::pair<int, int> PII;

inline char fgc() {
	static char buf[100000], * p1 = buf, * p2 = buf;
	return p1 == p2 && (p2 = (p1 = buf) + fread(buf, 1, 100000, stdin), p1 == p2)
		? EOF : *p1++;
}

inline LL readint() {
	LL res = 0, neg = 1; char c = fgc();
	for(; !isdigit(c); c = fgc()) if(c == '-') neg = -1;
	for(; isdigit(c); c = fgc()) res = res * 10 + c - '0';
	return res * neg;
}

inline char readsingle() {
	char c;
	while(!isgraph(c = fgc())) {}
	return c;
}

const int MAXN = 105, MAXM = 2005;
const int MO = 998244353;

int n, m;
LL a[MAXN][MAXM], s[MAXN], f[MAXN][MAXN << 1], g[MAXN][MAXN];

inline LL A(int in, int im) { // 让方程不那么长的简化
	return (s[in] - a[in][im] + MO) % MO;
}

int main() {
	n = readint(); m = readint();
	for (int i = 1; i <= n; i++) {
		for (int j = 1; j <= m; j++) {
			a[i][j] = readint();
			s[i] = (s[i] + a[i][j]) % MO;
		}
	}

	LL sum1 = 0, sum2 = 0;
	for (int im = 1; im <= m; im++) {
		memset(f, 0, sizeof(f));
		f[0][n] = 1;
		for (int in = 1; in <= n; in++) {
			for (int j = n - in; j <= n + in; j++) {
				f[in][j] = (f[in - 1][j] + f[in - 1][j - 1] * a[in][im] % MO + f[in - 1][j + 1] * A(in, im) % MO) % MO;
			}
		}
		for (int i = 1; i <= n; i++) {
			sum1 = (sum1 + f[n][n + i]) % MO;
		}
	}
	g[0][0] = 1;
	for (int i = 1; i <= n; i++) {
		for (int j = 0; j <= n; j++) {
			g[i][j] = g[i - 1][j];
			if (j > 0) g[i][j] = (g[i][j] + g[i - 1][j - 1] * s[i] % MO) % MO;
		}
	}
	for (int i = 1; i <= n; i++) {
		sum2 = (sum2 + g[n][i]) % MO;
 	}
	printf("%lld", (sum2 - sum1 + MO) % MO);
	return 0;
}