胡梦辰 朱滔 蒋青松 邹锐 吴桢 张晓玲 叶瑞 刘永,

滇池溶解氧浓度变化的氮磷循环响应模拟研究

胡梦辰1朱滔2蒋青松1邹锐3,4,5吴桢1张晓玲1叶瑞5刘永1,†

1.北京大学环境科学与工程学院, 国家环境保护河流全物质通量重点实验室, 北京 100871; 2.昆明滇池投资有限公司, 昆明 650034; 3.云南高原湖泊流域污染过程与管理重点实验室, 昆明 650034; 4.北京英特利为环境科技有限公司, 北京 100080; 5.南京智水环境科技有限公司, 南京 210012; †通信作者, E-mail: yongliu@pku.edu.cn

基于滇池外海的三维水动力‒水质模型, 模拟湖泊水动力及水质的动态过程, 识别湖泊水体中溶解氧特征, 并通过设置外部负荷削减情景, 探究溶解氧浓度对底泥内源释放以及氮磷循环的影响。结果表明: 1)由于藻类爆发与水体氧传输受阻的共同作用, 滇池底层水体 6—9 月间存在较严重的缺氧现象; 2)溶解氧调控底泥氮、磷营养物质的释放, 进而显著地影响底层水体中氮、磷浓度; 3)滇池底层缺氧状况的改善对外部负荷削减要求较高; 4)底泥的磷吸收‒释放平衡点溶解氧浓度约为 3.3mg/L, 滇池水体中氮浓度的变化更易受外源负荷作用影响。

水质模型; 溶解氧; 氮磷循环; 底泥释放; 滇池

影响湖泊水质的氮、磷营养盐主要来自外源、大气沉降和内源等[1–2], 在湖泊外部输入负荷不断削减的情况下, 底泥可能转化为湖泊水体中污染物的重要来源[3], 底泥参与的内源释放及营养物内循环持续地阻碍湖泊水质净化和生态恢复[3–4]。溶解氧(dissolved oxygen DO)对评估湖泊富营养化具有重要作用, 一方面, 作为重要的水质评价指标, 直接反映湖泊水体的污染程度; 另一方面, 作为水体组分, 广泛地影响各种生物地球化学循环过程[5–6]。同时, 在沉积物‒水体界面, 底层水(沉积物上覆水)的溶解氧作为重要的驱动因素, 对底泥内源释放和物质循环起着重要的调控作用[7–11]。因此, 湖泊富营养化相关研究中, 需要关注水体溶解氧浓度与氮磷循环之间的关系, 尤其是水体溶解氧浓度对湖泊底泥内源释放的响应规律。

研究湖泊溶解氧浓度对氮磷循环响应的方法主要包括监测数据分析、实验研究和湖泊模型模拟[5–17]。通过监测数据分析和实验研究, 在局部环境的机理阐释方面已有较多的发现。例如, Höhener等[9]发现, 溶解氧浓度是调控沉积物成为无机氮的汇或源的重要条件; Bertram[15]借助伊利湖中心区域1970—1986 年的长期监测数据, 发现水体磷酸盐含量与溶解氧浓度负相关; Wang 等[10]采集太湖水体和底泥样品, 通过实验识别厌氧环境对磷组分之间相互转化的影响机制。湖泊模型可实现对水体中各种反应和过程的定量刻画, 从而在改变边界条件的前提下, 更好地探索溶解氧浓度对氮磷循环的响应关系。例如, Bocaniov 等[16]开发了三维水动力学‒水生态模型 CAEDYM, 对不同强度磷入湖负荷情景下伊利湖中央盆地的缺氧现象进行评估; 邹锐等[17]基于三维水动力‒水质模型, 核算滇池氮磷通量的年变化, 揭示湖体内发生的重要转化过程。

基于以环境流体动力学模型(the environmental fluid dynamics code, EFDC)为内核的三维水动力‒水质模型[18], 本文以滇池外海为研究对象, 建立湖体底层溶解氧浓度对底泥释放及湖体水质的响应关系。通过不同外部条件下的模型模拟, 获取溶解氧、氮、磷浓度以及底泥释放通量等指标的年内变化过程, 分析溶解氧浓度的变化特征及其参与氮磷循环的规律, 以期为理解滇池的内循环过程提供参考。

1 研究对象和研究方法

滇池是我国重点治理的“三湖”之一, 位于云南省昆明市城区西南部, 是昆明市沿湖地区唯一的纳污水体。滇池流域面积为 2920km2, 以约占云南省0.75%的土地面积, 承载全省约 23%的 GDP 和 8%的人口。湖体北部有一道天然堤坝, 将滇池水域分割成南北两部分。南部称为外海, 是滇池的主体, 湖面面积为 298km2[1], 多年平均水资源量为 8.8 亿m3, 占滇池水资源总量的 90%以上。

本文的研究对象为滇池外海, 地处低纬度、高海拔地区, 属亚热带高原季风气候, 日温差较大, 冬、春季晴朗少雨, 夏、秋季湿润多雨, 5—10 月的降水量占全年总降水量的 70%～75%。滇池流域水资源缺乏, 生态用水不足, 内源累积性污染严重。滇池属宽浅型半封闭高原湖泊, 换水周期长, 自净能力差, 水生态系统退化严重, 易发生蓝藻水华。快速的城市扩张和工业发展, 对滇池造成高强度、长时期的污染, 底泥中积累大量污染物。

基于经过校验的滇池三维水动力‒水质模型[18], 可实现对湖泊中主要过程、水动力和水质指标的表征。进一步地, 结合模型中的生物和底泥等模块, 可实现对全湖任意网格和特定时间下的浓度数值输出。该模型采用离散化网格进行数值模拟, 需要考虑各个网格中水质组分的时空分异性, 各相关状态变量的存量通过对全湖所有网格进行积分获取[17]:

(,)=∭(,)ddd, (1)

其中,为组分的存量,为湖体中组分在时刻的浓度,,和为对应网格单元的空间特征坐标。通过对三维坐标求积分, 实现对全部湖泊单元网格的遍历, 从而输出时刻组分的全湖存量数据。

为了重点描述底泥作为内源在滇池中随底层溶解氧浓度变动而发生的变化, 对湖体底泥‒水界面的平面网格进行积分, 获取全湖底泥释放氮磷的累计通量[17]:

Ben(i)=∭b(i)ddd, (2)

其中,Ben(i)为累积底泥‒水体通量,b(i)为组分的底泥释放速率。首先通过计算获得累积通量, 然后按照运行步长对相邻累积通量做差, 即可得到对应时段单位步长内产生的过程通量, 进而计算各物质的通量释放速率。

设置模型的基础模拟情景为 2010 年, 模拟结果输出的时间间隔为 0.5 天。由于滇池的水体置换周期较长, 初始条件的影响在模拟过程中可以持续较长时间[18]。对某一确定的情景设置, 连续运行 5年, 以便消除初始条件的影响, 模型达到稳定状态后, 即得到可靠的输出结果。对于各项水质指标, 选取外海 8 个国控断面(图 1)的输出结果, 取平均值作为浓度数据; 对于各物质通量, 分别输出氮磷循环累积通量结果, 取相邻时间的差值作为该时段的通量数据。

2 结果与讨论

2.1 底层溶解氧浓度的年内变化

根据模型输出结果, 得到基准条件下滇池溶解氧浓度的分布情况(图 2)。整体而言, 滇池溶解氧浓度在垂直方向上存在较大差异, 将模拟输出分为底层(沉积物表层及其以上 1/6 深度的水体)和其他层(底层以上的其他水体)。在时间维度上, 由于富氧作用和滇池富营养化过程的共同影响, 滇池溶解氧浓度整体上呈现以年为单位的周期性变化, 冬、春季逐渐上升, 3 月开始, 随着温度升高而降低; 夏季藻类爆发使得底层氧消耗增大, 溶解氧浓度逐渐下降, 并在 5—9 月呈现明显的低位波动。在空间维度上, 滇池溶解氧浓度的垂向分布随深度变化差异显著, 突出体现在藻类爆发时期的溶解氧浓度随深度增加呈现明显下降的趋势。

图1 滇池湖体监测点位分布

湖泊缺氧(hypoxia)现象, 一般被界定为溶解氧浓度小于 2.0mg/L[16,19]。滇池底层水体 6—9 月间整体上处于低氧水平, 其溶解氧浓度在较长时间内低于 2.0mg/L, 说明滇池存在较严重的底层缺氧状况, 而其他层没有出现持续缺氧现象。这是由两方面原因造成的, 一方面, 藻类大量爆发后, 其残体的分解消耗大量溶解氧; 另一方面, 在夏季光照强烈时, 表层水体的温度随之升高, 但底层水难以被加热, 此时滇池底层平均水温为 16.9℃, 而其他层平均水温为 23.6℃, 作为大型高原浅水湖泊, 以较小的水深变化呈现较大的水温差异, 垂向上的密度差影响氧的完全交换[20–21]。

2.2 氮磷浓度及通量波动

通过分析湖泊底层水体中磷酸盐(PO4)、氨氮(NH3-N)和硝氮(NO3-N)与溶解氧浓度的关系(图 3), 并核算底泥内源释放的氮磷物质通量(图 4, 通量为负值时, 表示底泥从水体中吸收氮磷营养物质), 研究溶解氧浓度的变化及其对底泥源汇功能的影响。

在大部分时期, 磷酸盐浓度与溶解氧浓度负相关, 但存在异常的月份区间, 这一现象由磷释放机制和滇池富营养化过程共同造成。从图 3 看出, 6 月—次年 2 月, 溶解氧浓度与磷酸盐浓度明显负相关, 底泥中大部分磷以矿物结合态存在, 随着溶解氧浓度提高, 底层氧化氛围加强, 使得磷酸盐也被铁、锰的高价态化合物吸附并沉降, 水体中浓度降低[6,22]。在低溶解氧浓度下, 这一可逆过程的方向反转, 底泥中的磷酸盐释放到水体中。3—5 月, 随着温度升高, 湖泊中藻类生长加快, 底层溶解氧被快速消耗, 同时藻类生物量的快速增加使得藻类对磷酸盐的吸收超过溶解氧驱动, 成为这一阶段磷酸盐变化的主要原因, 此时溶解氧浓度与磷酸盐浓度均呈现下降趋势。

氨氮和硝氮浓度在全年时间尺度上呈现一致的趋势(图 3): 氨氮浓度与溶解氧浓度负相关, 硝氮浓度与溶解氧浓度正相关。结合底泥释放含氮物质过程可知, 溶解氧作为反应物参与硝化作用, 促进氨氮等含氮物质转换为硝氮; 同时, 高浓度的溶解氧抑制湖泊与底泥中反硝化过程, 减少水体中硝氮的消耗[4–5]。如图 4 所示, 底泥氨氮释放通量与溶解氧浓度负相关, 而硝氮释放通量与溶解氧浓度呈现较好的正相关关系, 尤其在 6—9 月的低溶解氧水平下, 底泥大量释放氨氮, 而硝氮释放相对减缓。

图2 滇池溶解氧浓度的年内变化

图3 滇池底层水体中氮磷营养物质浓度的年内变化

图4 滇池氮磷营养物质底泥释放通量的年内变化趋势

2.3 不同负荷削减情景下溶解氧浓度的变化

利用三维水动力‒水质模型[18], 模拟滇池在边界条件变化下的状态(尤其是无法进行实际观测的情景), 探究需要极端变化才能体现的特定系统响应。外源营养物质负荷削减是富营养化控制的重要手段[10,23], 能够改善湖泊的缺氧现象。本文选择不同的外源负荷削减情景, 基于模型得到不同削减强度下滇池水体溶解氧浓度的垂向分布情况。图 5 显示, 低、中强度的外部负荷削减对改善滇池水体底层缺氧现象的效果微弱, 即使在 50%的削减强度下, 夏季滇池底层缺氧现象仍较严重, 削减强度达到90%时才有较明显的改善; 在藻类爆发期(5—10月), 外源负荷削减导致藻类生物量下降, 随削减强度增加, 表层水体溶解氧浓度出现较明显的下降, 使得向底层进行溶解氧传质的驱动力减弱, 反而不利于底层缺氧现象的改善; 外源负荷削减明显提高旱季(11 月—次年 4 月)湖泊水体的溶解氧浓度, 但旱季水体处于良好的溶解氧水平, 因此负荷削减的贡献不明显。

2.4 溶解氧浓度对氮磷通量变化的响应

从前面的讨论可知, 在底泥与底层水的界面, 各种形态的氮磷物质通量差异显著。为探究其成因, 选取基准、50%负荷削减和 90%负荷削减 3 种典型情景展开分析。图 6～8 以水体温度为色阶, 对比 3 种情景下底泥释放过程与溶解氧浓度的关系。

图 6 中, 外源负荷削减后, 数据点的分布更加集中, 基准情景下低、中、高温度的 3 个清晰边界转变为两个, 90%削减情景下, 温度的影响更加模糊。由于外部负荷的削减, 温度对通量‒溶解氧响应的影响大大减小, 散点包络围成的区域明显缩小, 使得溶解氧浓度成为底泥磷释放的主导因素。上述现象符合同类研究[10,22]揭示的一般性规律: 底泥磷酸盐释放速率随着溶解氧浓度的升高而降低, 并伴随迁移方向的改变, 从向水体释放转变为从水体中吸收, 达到吸收阈值后保持稳定。

图5 滇池溶解氧浓度的垂向分布

图6 滇池底层水体溶解氧‒磷酸盐浓度的散点分布

图7 滇池底层水体溶解氧‒氨氮浓度的散点分布

图8 滇池底层水体溶解氧‒硝氮浓度的散点分布

由于滇池水动力系统的特性, 溶解氧参与底泥磷酸盐释放过程受温度影响显著。随着外源负荷削减强度增大, 温度在响应关系中的影响被压缩, 使得溶解氧成为底泥磷酸盐释放的主导。这一规律由滇池的地理气候条件、湖泊水质状况和底泥组成决定。选取对数形式回归方法进行拟合, 释放速率()与溶解氧浓度()的拟合方程为

= ‒0.957 ln()+1.1331 (2=0.69, sig.<0.01)。 (3)

滇池底泥磷酸盐的吸收与释放过程在底层溶解氧浓度约为 3.3mg/L 时达到动态平衡。由于滇池底层水体溶解氧浓度年内变化较大, 促使底泥出现源汇功能的动态转变: 6—9 月整体上处于低氧水平, 底泥作为内源, 向湖泊水体释放磷酸盐; 其他时间溶解氧浓度正常, 底泥作为汇, 从水体中吸收磷酸盐。滇池底泥源‒汇转变的规律导致外部负荷削减的效益在一定程度上被抵消, 增大了湖泊治理及水质达标的难度。

溶解氧浓度对底泥氮释放的影响随负荷削减强度的变化更加剧烈, 说明氮释放更容易受外部条件影响。

对于氨氮通量, 如图 7 所示, 一方面, 低溶解氧浓度可直接抑制底泥‒水体界面的硝化作用, 减少氨氮的消耗, 另一方面, 低氧环境下水体的还原性增强, 厌氧菌分解有机物产生氨氮, 增加氨氮的释放。温度更显著地影响厌氧菌的分解活动, 因此在同一溶解氧浓度下, 夏季高温环境中氨氮的释放强度相对更高。外源负荷降低的影响则体现在底泥氨氮释放过程中压缩了温度的影响, 随着负荷削减强度增大, 图 7 中散点包络区域缩小, 氨氮释放速率与底层溶解氧浓度负相关关系的线性增强, 低溶解氧浓度下最大氨氮释放强度减小。基准情景下, 夏季藻类爆发, 底层溶解氧下降, 底泥氨氮释放速率迅速上升, 向水体大量输出氮营养物质, 使得湖体水质进一步恶化。但是, 在高负荷削减强度下, 因藻类爆发导致缺氧环境, 底泥的氨氮释放强度也相对减弱, 与外部负荷削减协同促进水体氨氮浓度的降低。

对于硝氮通量, 如图 8 所示, 随溶解氧浓度增加, 底泥释放速率先升后降, 并且在两个阶段都呈现较好的线性关系。温度的影响体现在折线拐点上, 低温区间拐点出现较早, 高温区间拐点出现较晚。一方面, 溶解氧作为反应物参与硝化反应, 促进底泥释放硝氮; 另一方面, 底泥向水体的硝氮释放过程受表面传质系数影响, 而表面传质系数与溶解氧浓度成反比[24], 使得较高浓度的溶解氧反而抑制底泥硝氮的释放。温度同样影响表面传质系数, 高温下表面传质系数增大, 因此在高温区间, 需要更高的溶解氧浓度才能产生抑制作用。

在底泥硝氮释放过程中, 外源负荷的削减抑制硝氮释放强度, 却增强温度对溶解氧‒硝氮释放的影响。一方面, 硝氮释放年变化波动更加剧烈; 另一方面, 温度对表面传质系数的影响被放大, 使得在低温与高温两个区间, 底泥硝氮释放拐点对应的溶解氧浓度差距增大。90%负荷削减情景下, 低温区间拐点出现更早, 对应溶解氧浓度降低; 高温区间拐点出现更晚, 对应溶解氧浓度升高。

3 结论

本研究选择滇池外海作为研究对象, 基于滇池三维水动力‒水质模型, 模拟湖泊水体中各水质组分浓度及氮磷循环过程中通量的动态变化, 识别溶解氧浓度与氮磷浓度的关系, 主要结论如下。

1)虽然滇池属于浅水湖泊, 但也存在较显著的底层缺氧现象, 集中出现在 6 月至 9 月中旬, 是藻类爆发与水体氧传输受阻共同造成的。

2)溶解氧广泛参与湖泊的氮磷营养物质循环过程, 通过调控底泥氮磷释放过程而显著地影响底层水体中氮磷物质浓度。对于氨氮和磷酸盐的底泥释放, 溶解氧浓度的年内变化使得滇池底泥发生源‒汇功能的转变。

3)入湖负荷削减对湖泊溶解氧指标有改善作用, 但中低强度的负荷削减改善作用有限。负荷高度削减情景下, 溶解氧对氮磷释放通量的驱动作用更显著。

4)滇池底泥的磷吸收‒释放平衡点溶解氧浓度约为 3.3mg/L, 说明氮浓度更易受外源负荷作用的影响。

[1] 贺克雕. 滇池水质状况综合评价及变化趋势分析. 人民长江, 2012, 43(12): 37–41

[2] 刘永, 邹锐, 郭怀成, 等. 智能流域管理研究. 北京: 科学出版社, 2012

[3] Marsden M W. Lake restoration by reducing external phosphorus loading: the influence of sediment phos-phorus release. Freshwater Biology, 2010, 21(2): 139–162

[4] 张路, 范成新, 王建军, 等. 太湖水土界面氮磷交换通量的时空差异. 环境科学, 2006, 27(8): 1537–1543

[5] 殷燕, 吴志旭, 刘明亮, 等. 千岛湖溶解氧的动态分布特征及其影响因素分析. 环境科学, 2014, 35 (7): 2539–2546

[6] Robertson D M, Imberger J. Lake number, a quanti-tative indicator of mixing used to estimate changes in dissolved oxygen. International Review of Hydrobio-logy, 2010, 79(2): 159–176

[7] Nowlin W H, Evarts J L, Vanni M J. Release rates and potential fates of nitrogen and phosphorus from sedi-ments in a eutrophic reservoir. Freshwater Biology, 2005, 50(2): 301–322

[8] Foley B, Jones I D, Maberly S C, et al. Long-term changes in oxygen depletion in a small temperate lake: effects of climate change and eutrophication. Freshwater Biology, 2012, 57(2): 278–289

[9] Höhener P, Gächter R. Nitrogen cycling across the sediment-water interface in an eutrophic, artificially oxygenated lake. Aquatic Sciences, 1994, 56(2): 115–132

[10] Wang S, Jin X, Bu Q, et al. Effects of dissolved oxygen supply level on phosphorus release from lake sediments. Colloids and Surfaces A: Physicochemical and Engineering Aspects, 2008, 316(1): 245–252

[11] Hupfer M, Lewandowski J. Oxygen controls the phosphorus release from lake sediments —a long-lasting paradigm in limnology. International Review of Hydrobiology, 2010, 93(4/5): 415–432

[12] Zhou Y, Obenour D R, Scavia D, et al. Spatial and temporal trends in Lake Erie Hypoxia, 1987–2007. Environmental Science & Technology, 2013, 47(2): 899–905

[13] Berkheiser V E, Street J J, Rao P S C, et al. Partitioning of inorganic orthophosphate in soil-water systems. Critical Reviews in Environmental Science and Technology, 1980, 10(3): 179–224

[14] Moore P A, Reddy K R, Graetz D A. Nutrient trans-formations in sediments as influenced by oxygen sup-ply. Journal of Environmental Quality, 1992, 21(3): 387–393

[15] Bertram P E. Total phosphorus and dissolved oxygen trends in the central basin of Lake Erie, 1970–1991. Journal of Great Lakes Research, 1993, 19(2): 224–236

[16] Bocaniov S A, Leon L F, Rao Y R, et al. Simulating the effect of nutrient reduction on hypoxia in a large lake (Lake Erie, USA-Canada) with a three-dimen-sional lake model. Journal of Great Lakes Research, 2016, 42(6): 1228–1240

[17] 邹锐, 吴桢, 赵磊, 等. 湖泊营养盐通量平衡的三维数值模拟. 湖泊科学, 2017, 29(4): 819–826

[18] Wang Z, Zou R, Zhu X, et al. Predicting lake water quality responses to load reduction: a three-dimen-sional modeling approach for total maximum daily load. International Journal of Environmental Science & Technology, 2014, 11(2): 423–436

[19] Nürnberg G K. Quantifying anoxia in lakes. Limno-logy & Oceanography, 1995, 40(6): 1100–1111

[20] Li W C, Yang X D, Yin Y, et al. Thermal stratification in lake Zige Tangco, central Tibetan Plateau. Wuhan University Journal of Natural Sciences, 2005, 10(4): 689–693

[21] Snortheim C A, Hanson P C, Mcmahon K D, et al. Meteorological drivers of hypolimnetic anoxia in a eutrophic, north temperate lake. Ecological Model-ling, 2017, 343: 39–53

[22] 梁文, 王泽, 焦增祥, 等. 内源磷的释放作用及影响因素研究进展. 四川环境, 2012, 31(5): 105–109

[23] Liu Y, Wang Y, Sheng H, et al. Quantitative evalu-ation of lake eutrophication responses under alterna-tive water diversion scenarios: a water quality mo-deling based statistical analysis approach. Science of the Total Environment, 2014, 468/469(2): 219–227

[24] Ji Z G. Hydrodynamics and water quality: modeling rivers, lakes, and estuaries. Hoboken: Wiley, 2008

Simulation Study on Nitrogen and Phosphorus Reycling Response of Changing Dissolved Oxygen Concentration in Lake Dianchi

HU Mengchen1, ZHU Tao2, JIANG Qingsong1, ZOU Rui3,4,5, WU Zhen1, ZHANG Xiaoling1, YE Rui5, LIU Yong1,†

1. College of Environmental Science and Engineering, State Environmental Protection Key Laboratory of All Materials Flux in Rivers, Beijing 100871; 2. Kunming Dianchi Investment Co., Ltd, Kunming 650034; 3. Yunnan Key Laboratory of Pollution Process and Management of Plateau Lake-Watershed, Kunming 650034; 4. Beijing Inteliway Co., Ltd, Beijing 100080; 5. Nanjing Innotwater Environmental Technology Co., Ltd, Nanjing 210012; † Corresponding author, E-mail: yongliu@pku.edu.cn

A three-dimensional hydrodynamic-water quality model of Lake Dianchi was used to simulate the dynamic process of lake hydrodynamics and water quality, identify the characteristics of dissolved oxygen (Do) in the lake water. Combined with the setting of the external load reduction scenario, the effects of dissolved oxygen on the endogenous release of sediment and the cycle of nitrogen and phosphorus were explored. The results are conckyded. 1) The severe hypoxia in Lake Dianchi from June to September is caused by a combination of algae outbreaks and obstructed oxygen transport in the water. 2) Dissolved oxygen regulates the release of nitrogen and phosphorus nutrients in the sediment, and then significantly affects the concentration of nitrogen and phosphorus in the bottom water. 3) The improvement of the hypoxia condition at the bottom of the lake has high requirements on the reduction of external load. 4) The dissolved oxygen concentration at the phosphorus absorption-release equilibrium point of the sediment is about 3.3 mg/L, and the change in nitrogen concentration in Lake Dianchi is more susceptible to external load.

water quality model; dissolved oxygen (DO); nitrogen and phosphorus cycle; sediment release; Lake Dianchi

10.13209/j.0479-8023.2021.024

“十三五”规划项目“滇池流域水污染控制工程评估及精准治污决策系统研究”(JKM2019070419)、国家自然科学基金(51779002)和水利部珠江河口动力学及伴生过程调控重点实验室开放研究基金(No.[2017]KJ11)资助

2020–04–05;

2020–07–29