复杂物质循环过程模拟方法与平台实现——以畜禽养殖系统磷循环为例

删除或更新信息，请邮件至freekaoyan#163.com(#换成@)

本站小编 Free考研考试/2021-12-29

盛虎^,¹, 刘欣², 芦昕雨², 廖宇星², 袁增伟^,¹^,²

1.南京大学（溧水）生态环境研究院,南京 211200

2.污染控制与资源化研究国家重点实验室,南京大学环境学院,南京 210023

Complex material cycling process simulation method and platform: Taking the phosphorus cycling of livestock and poultry breeding system as an example

SHENG Hu^,¹, LIU Xin², LU Xinyu², LIAO Yuxing², YUAN Zengwei^,¹^,²

1. Lishui Institute of Ecological Environment, Nanjing University, Nanjing 211200, China

2. State Key Laboratory of Pollution Control and Resource Reuse, School of the Environment, Nanjing University, Nanjing 210023, China

通讯作者: 袁增伟,男,河南平顶山人,教授,博士生导师,研究方向为物质循环及其环境效应。E-mail: yuanzw@nju.edu.cn

收稿日期:2020-11-4修回日期:2021-03-9网络出版日期:2021-03-25

基金资助:

国家重点研发计划项目.2018YFF0215503
国家****科学基金项目.41925004

Received:2020-11-4Revised:2021-03-9Online:2021-03-25
作者简介 About authors
盛虎,男,湖北孝感人,博士,副研究员,研究方向为磷循环过程模拟及环境大数据。E-mail: shenghu@nju.edu.cn

摘要
物质循环是人类获取资源的重要方式,尤其是对于磷这种地球生物生存不可或缺、不可替代的资源,理解其循环过程对资源可持续利用乃至人类永续发展具有重要意义。目前关于物质循环的研究多采用物质流分析方法,该方法是一种以物质流动过程质量守恒为基础、以系统分析为手段的定量核算方法,需要将系统物质流动路径与核算模型融为一体。为了解决传统物质流分析方法中数据与模型无法分离从而不适应与复杂物质循环过程模拟的问题,本文提出了复杂物质循环过程模拟概念模型与定量计算公式,建立了基于数据库、模型库、结果库、分析库4库分离的复杂物质循环过程建模方法,并基于EnVirLab环境虚拟仿真实验平台,以中国1600—2012年畜禽养殖系统磷循环过程为例开展方法验证研究,验证了基于4“库”分离的复杂物质循环过程模拟分析方法的可行性。本文对于建立标准化的物质流分析方法及开发复杂物质循环过程模拟软件,进而精确量化物质循环路径,对实现物质循环过程模拟从定量格局描述性分析向定量机理解释性分析转变具有重要意义。
关键词： 物质流分析;磷循环;复杂网络;虚拟仿真;系统动力学;软件平台

Abstract

Material cycling is an important way for human beings to obtain resources, especially for phosphorus, an indispensable and irreplaceable resource for life on the Earth, and it is of great significance to understand phosphorus cycling process for the sustainable use of resources and sustainable development of human society. At present, the research on material cycle uses the material flow analysis method, which is a quantitative accounting method based on the quality conservation of the material flow process and by means of systematic analysis, integrating the system flow path with the accounting model. In order to solve the problem that the data and model in the traditional material flow analysis method cannot be separated, which is not suitable for the simulation of complex material cycling processes, this study proposed the conceptual model and quantitative calculation formula of complex material cycle process simulation, and established the complex material cycling process modeling method based on database, model library, result library, and analysis library separation. Based on the environment virtual laboratory (EnVirLab) platform and taking the phosphorus cycling process of livestock and poultry breeding system in China from 1600 to 2012 as an example, a method validation study was carried out, and the feasibility of complex material cycling process simulation and analysis based on four library separation was verified. This study is of significance for establishing standardized material flow analysis methods and developing complex material cycling process simulation software, quantifying the material cycling path accurately, and realizing the transformation of material cycling process simulation from descriptive quantitative pattern analysis to explanatory quantitative mechanism analysis.

Keywords：material flow analysis;phosphorus cycling;complex network;virtual simulation;system dynamics;software platform

PDF (11812KB)元数据多维度评价相关文章导出EndNote|Ris|Bibtex 收藏本文
本文引用格式
盛虎, 刘欣, 芦昕雨, 廖宇星, 袁增伟. 复杂物质循环过程模拟方法与平台实现——以畜禽养殖系统磷循环为例. 资源科学[J], 2021, 43(3): 465-476 doi:10.18402/resci.2021.03.04
SHENG Hu, LIU Xin, LU Xinyu, LIAO Yuxing, YUAN Zengwei. Complex material cycling process simulation method and platform: Taking the phosphorus cycling of livestock and poultry breeding system as an example. RESOURCES SCIENCE[J], 2021, 43(3): 465-476 doi:10.18402/resci.2021.03.04

1 引言

随着人类社会经济的飞速发展,大量的物质资源从自然界进入人类社会,并在人类社会生产、分配、交换和消费各个过程之间迁移与转化^[1,2]。其中,磷是支撑地球上生命系统持续发展的一种不可替代的资源,其循环过程对于人类的生存和发展具有重要意义^[3,4]。然而,人口增长和生活水平持续提高导致粮食需求量剧增,人类开始大量开采磷矿资源并严重改变了自然界磷循环路径和磷素流动强度,进而引发了磷矿资源枯竭^[5,6]和水体富营养化^[7]等一系列问题。在行星边界理论（Planetary Boundary, PB）中,磷已经成为影响人类社会持续发展的重要因素,并在很大程度上决定着地球系统是否能够保持在弹性范围和适应状态^[8],事实上其影响正随着时间不断被突破,目前已经进入到高风险区域^[9]。磷循环过程加剧会导致至少3个方面的后果^[10,11]：第一,由于磷矿石大量开采造成磷矿资源化枯竭风险,而这必然影响到磷肥供应并进而带来全球粮食产量下降风险,有可能加剧全球粮食危机;第二,大量的磷矿石被开采加工成磷肥以支持农业生产,但磷资源的利用效率不高,导致大量的磷流失到环境中,尤其是江河湖泊等水环境,造成水体磷浓度超标;第三,水环境中的磷经过一系列生物化学过程后,引起水中藻类大量繁殖从而导致水体发生富营养化,最终威胁水生态系统生物多样性。可见,无论是从资源角度,还是环境角度,磷资源不合理开发利用对于人类以及生态环境的影响都是十分严重的,甚至是不可逆的。那么,研究如何有效利用类似于磷这样的物质资源并减少其向自然界中排放具有重要意义^[12,13,14]。然而,在日益复杂的人类活动背景下,物质循环过程的复杂度也随之增加^[15],具体表现为：①人类活动的加剧使得物质循环过程的变化也逐渐加剧,从而导致研究的时间精度需要不断细化,研究的时间长度需要不断延伸;②由于分工协作与集约化发展使得人类活动的空间异质性增强,从而导致研究的空间分辨率也逐渐精细化;③随着人类活动的细化及对其研究分析的不断深入,物质循环过程中涉及的要素与过程及其相互关系也趋于复杂。这些因素导致物质循环过程模拟的复杂度也不断提升,一方面阻碍了对物质循环过程的理解,另一方面也加大了模拟模型的计算量,从而导致物质循环过程模拟缺乏有效的建模与求解方法。为此,构建一套复杂物质循环过程模拟方法并开发模拟模型求解计算平台十分必要。

关于物质循环过程的模拟,通常采用的是物质流分析方法^[16,17]（Substance Flow Analysis, SFA）。该方法基于物质守恒原理,通过系统分析的方法确定系统边界内物质循环的各个过程（系统要素）,以及各个过程之间的物质流（要素联系）,采用定量核算的方法（如系数法）来确定物质流的大小^[18,19,20]。具体包括两个主要步骤：①定性分析,包括系统边界划定、系统过程确定、过程之间关系分析等;②定量分析,包括核算方法确定、数据采集、模型计算和结果分析等步骤。对于定性分析,通常采用Microsoft Visio、Microsoft PowerPoint、e!Sankey等流图绘制工具来构建概念模型^[21];而对于定量分析,则主要采用Spreadsheet的计算方法（比如Microsoft Excel等）进行计算和作图分析。这种分析方法将物质循环的分析过程割裂成两个独立的部分,容易造成分析思路不连贯以及概念模型与核算模型不一致的问题,而且随着研究系统复杂度的增加,出错的可能性也会逐渐增大。为了解决这个问题,维也纳工业大学水质、资源与废物管理研究所研究团队开发了STAN^[22]（subSTance flow ANalysis的简称）建模软件平台。该软件平台实现了流图绘制与模型核算结果的一体化,在可视化方面可以通过线条的粗细来表征流的大小,并且支持Microsoft Excel数据导入导出。STAN自开发以来被广泛地应用于物质流（如磷流）分析中^[23,24,25]。然而STAN也并非是一个真正的核算平台,而仅仅只是一种可视化与数据校正（过程的平衡分析）工具。随着对物质循环过程理解和研究的进一步深入,物质循环过程的定量刻画逐渐变得更加复杂,模拟的时间长度和空间分辨率要求也越来越高,数据结果的可靠性和不确定性也成为建模分析的重要评估指标。在这个背景下,采用传统的基于Spreadsheet的计算方法显然已经不能满足当前的研究需求。针对这一问题,在研究中国1600—2012年磷循环过程时,Liu等^[4]提出了一种新的物质流核算技术,该技术将复杂的概念模型以矩阵形式加以刻画,同时将核算模型从复杂的Spreadsheet中提取出来,实现模型与数据的分离。按照这种思路开发了R软件包sfc^[26]（Substance Flow Computation）,成功实现了长时间序列复杂磷循环过程的定量建模分析和不确定性评估。然而,该核算技术脱离了流图的概念,直接通过抽象的数学公式来刻画模型,需要建模者具备一定抽象分析能力,难以推广使用。

由此可见,目前已有的分析方法与分析工具尚不能很好地解决复杂物质循环过程模拟的问题,主要表现在可视化的建模工具难以支撑复杂模型的建模与求解,而复杂模型的建模分析工具又缺乏可视化的建模手段。出现这些问题的原因是目前对于复杂物质循环过程的认识不够,对于定量表达与求解复杂物质循环过程模拟模型缺乏方法与工具支撑。为此,基于前期的研究工作,本文试图提出一种复杂物质循环过程模拟模型及其平台实现方法,并通过中国畜禽养殖系统磷循环模拟开展实证研究,以期回答3个方面的问题：①什么是复杂物质循环过程?②如何定量表达与求解复杂物质循环过程?③如何验证这种方法的可行性?

2 复杂物质循环过程

本文所说的“复杂物质循环过程”,是指在研究物质循环过程时,由于时间、空间及过程3个维度颗粒度及范围提升导致系统复杂度急剧增加的现象。这种现象可以通过图1所示的复杂物质循环过程概念框架加以理解。图1包括3个维度,其中横坐标表示的是空间维度,纵坐标表示的是时间维度,方格表示的是过程维度。每个区域在每个时间点上都发生一系列人类活动过程,如果把这些过程分成系统1、系统2、…、系统n,这些过程之间存在物质的相互流动,同时这些过程会从环境中获取相关物质资源,也会将含有物质的废物排放到环境中。此外,这些过程通过贸易方式会与其他区域的活动过程发生物质的相互流动,每个过程可能会产生物质的存量,这些存量会随着时间的变化从一个时间点传递到另一个时间点。物质在迁移转化的过程中,满足物质守恒原理。如果将整个复杂物质循环过程抽象成一个有向图,那么这个有向图就是一个包含时间和空间信息的复杂网络。

图1

新窗口打开|下载原图ZIP|生成PPT
图1复杂物质循环过程概念框架

Figure 1Conceptual framework of complex material cycling process

对于复杂物质循环过程,需要建立一种简便易行的分析方法。如果图1中每个方格里面涉及的过程都不相同,而且这些过程之间的连接也复杂多样,那么几乎不可能把这样一个复杂问题简化。事实上,图中每个方格在结构上是十分相似的。其原因是因为同一时间上不同区域的人类活动相似,或者同一区域在一定时间范围内不同时间点上的人类活动由于“惯性”或“路径依赖”的缘故也相似。因此,可以先忽略掉“时间”与“空间”维度,而只考虑“过程”这个维度。而对于“过程”维度,又可以将其拆分成不同的“系统”。按照这种思路,复杂物质循环过程的分析问题就能拆解成简单的系统过程流分析问题：对每个“系统”建立核算模型,并通过“系统”与“系统”之间的物质流,将这个核算模型串成“全过程”核算模型,最后将以上分析过程应用在不同的“时间”维度和“空间”维度。至于串接方式,在时间维度上,通过系统存量的变化来串接;在空间维度上,通过区域之间的贸易来串接。

由此可见,需要建立一种方法实现任何一个区域内部、任何一个系统长时间序列的物质流动态模拟,并通过移植和复制的方法对所有区域、所有系统的物质流进行动态刻画。

3 模型构建与求解方法

3.1 数学模型构建

由前面的分析可知,解决复杂物质循环过程的分析问题,关键就是构建单系统物质流动态建模的方法。物质流分析方法是一个很好的切入点,因为它是定量刻画物质循环过程的一种基本方法,其原理就是通过系统分析来将复杂系统划分成不同的过程,定量刻画过程与过程之间的物质流动关系。对于每一个过程,其输入流与输出流满足物质守恒原理。然而,物质流分析方法仅仅提供了一种分析框架,并没有提出该分析框架的具体实现方法。为此,本文提出了一种针对单系统细化物质流分析框架的方法（图2）。

图2

新窗口打开|下载原图ZIP|生成PPT
图2单系统物质流动态模拟方法原理

Figure 2Principle of single system material flow dynamic simulation method

图2中①号框的内容就是单系统物质流动态模拟概念模型。①号框内主要研究的是“中间过程”,它的输入流来源于“上游过程”,输出流来源于“下游过程”。在实际建模过程中,“上游过程”和“下游过程”可能不止一个,此处为了简化仅仅给出了一个。记输入流为

I (t)

,输出流为

O (t)

,这里

t

表示的是第

t

年。由于输入流和输出流可能不止一条（图2中的输入流和输出流都是2条）,那么每一条输入流和输出流可以分别命名为

I_{i} (t)

和

O_{j} (t)

,其中

I_{i} (t)

表示第

i

项（

i = 1, ?, n

）输入在第

t

时刻的输入流;

O_{j} (t)

表示第

j

项（

j = 1, ?, m

）输出在第

t

时刻的输出流。对于图2中①号框,有如下定量关系：

（1）

I (t) = \overset{n}{\sum_{i = 1}} I_{i} (t)

（2）

O (t) = \overset{m}{\sum_{j = 1}} O_{j} (t)

对于“中间过程”,同时还需要满足物质守恒的条件,即：

（3）

Δ S (t) = I (t) - O (t)

式中：

Δ S (t)

表示系统存量在第

t

时刻相比于在第

t - 1

时刻的增加量。在稳态条件下,

Δ S (t) = 0

,即系统输入输出平衡。

公式（1）-（3）即为物质流分析方法的定量表达,该表达形式需要进一步细化,通过图2中②号框的内容来表达。②号框内主要研究的是每一条流的计算方法,在此仅仅以“输入流1”为例开展说明：“输入流1”是由“变量1”和“参数1”计算得出（图2中通过虚线箭头表示这种计算关系,箭头方向为自变量指向因变量）,“变量1”可以是一个变量,也可以是一系列变量,“参数1”可以是一个参数,也可以是一系列参数。变量与参数的差别在于,变量是系统内部随时间而变化的量,而参数则是系统的固有属性,不随时间发生变化。对于“中间过程”而言,如果假设输入流的变量为

A_{i} (t)

,参数为

α_{i}

;输出流的变量为

B_{j} (t)

,参数为

β_{j}

,计算函数分别为

f_{I}^{i}

和

f_{O}^{j}

,那么有如下计算公式：

（4）

I_{i} (t) = f_{I}^{i} (A_{i} (t), α_{i})

（5）

O_{j} (t) = f_{O}^{j} (B_{j} (t), β_{j})

公式（4）-（5）为“中间过程”的具体计算公式。“中间过程”只是所有过程中的一个特例,考虑到记号

I

和

O

分别表示输入和输出,而对于“中间过程”的“输入流”同时也是“上游过程”的输出流,记任意中间过程的编号为

P_{k}

,其上游过程的编号为

P_{k - 1}

,其下游过程为

P_{k + 1}

,那么记号

I

可以表示为

P_{k - 1} \to P_{k}

,记号

O

可以表示为

P_{k} \to P_{k + 1}

。另外,

i

和

j

统一规约为s,

A

和

B

统一规约为

X

α

和

β

统一规约为

θ

。那么公式（4）-（5）可以重新记为：

（6）

F_{s}^{P_{k - 1} \to P_{k}} (t) = f_{P_{k - 1} \to P_{k}}^{s} (X_{s} (t), θ_{i})

式中：

F_{s}^{P_{k - 1} \to P_{k}} (t)

即为

I_{s} (t)

;

F_{s}^{P_{k} \to P_{k + 1}} (t)

即为

O_{s} (t)

。通过公式（6）,可以计算出所有过程的输入输出流。

根据以上分析,不难发现图2中方框①的内容是对物质循环过程的结构刻画,揭示了物质循环包含哪些过程,这些过程有哪些输入流和输出流,通过公式（1）-（3）进行定量表达;方框②的内容是对物质循环过程的流量刻画,揭示了每一条流与哪些变量和参数相关,是如何计算出来的,通过公式（4）-（5）或者公式（6）进行定量计算。

3.2 数值求解方法

数值求解复杂物质循环过程模拟模型是一项系统工程,为此本文提出了如图3所示的架构图来实现这一过程。从图3中可以看出,复杂物质循环过程模拟需要形成4个“库”,即数据库、模型库、结果库和分析库。之所以需要这样划分,是因为传统的计算方法往往把这些“库”结合在一起,从而增加了计算和分析的复杂度,不适合用于复杂计算。本文数据库的构建包括变量和参数两个部分,即公式（6）中的

X_{s} (t)

和

θ_{s}

。模型库包括概念模型和数学模型两个部分,即公式（6）中的

f_{P_{k - 1} \to P_{k}}^{s}

（其中概念模型由

P_{k - 1} \to P_{k}

来表征）。结果库包括图形和表格两个部分,即公式（6）中的

F_{s}^{P_{k - 1} \to P_{k}} (t)

。分析库主要开展不确定性分析和参数敏感性分析,通过Monte Carlo随机抽取参数

θ_{i}

然后代入模型

f_{P_{k - 1} \to P_{k}}^{s}

中,计算出结果

F_{s}^{P_{k - 1} \to P_{k}} (t)

,从而评估参数

θ_{s}

对结果

F_{s}^{P_{k - 1} \to P_{k}} (t)

的影响。

图3

新窗口打开|下载原图ZIP|生成PPT
图3复杂物质循环过程模拟架构图

Figure 3Framework of complex material cycling process simulation

值得注意的一个问题是：如何实现复杂物质循环过程模拟?对于这个问题的回答,首先需要对公式（3）作进一步分析。通常物质流分析是以年为时间尺度的,那么对于存量的变化,可以写成如下形式：

（7）

Δ S (t) = S (t) - S (t - 1) = (I (t) - O (t)) Δ t = 1

如果把时间尺度变成

Δ t

,那么公式（7）则可以变为

（8）

Δ S (t) = S (t) - S (t - Δ t) = (I (t) - O (t)) Δ t

即

（9）

\frac{Δ S (t)}{Δ t} = \frac{S (t) - S (t - Δ t)}{Δ t} = I (t) - O (t)

当

Δ t \to 0

时,公式（9）可以变为

（10）

\frac{dS (t)}{dt} = I (t) - O (t)

由此可见,物质流分析方法可以通过常微分方程组来进行表达。动态物质流分析就是按照这个思路开展相关研究的^[27]。本质上讲,长时序的物质流分析都是动态物质流分析。在求解微分方程组时,把时间步长设置为1（年）,就可以得到逐年的计算结果。

目前,求解常微分方程组最方便的工具就是用于求解系统动力学的分析工具。从方法的应用上看,动态物质流分析与系统动力学分析是两种不同的分析方法,解决的问题不一样,分析的思路也不一样,但是抛掉学科背景,其数学本质都是求解微分方程组,因此在求解上是可以借用的。常用的系统动力学分析软件,包括商业版软件（如Stella、Vensim、Powersim等）和开源免费版软件（如NetLogo、Simantics、Insight Maker等）。本文为了验证复杂物质循环过程模拟的可行性,在Insight Maker源代码的基础上开发了环境虚拟仿真实验平台EnVirLab（通过https://envirlab.netlify.app/网址即可以访问使用,可以在不保存数据和模型的条件下在线运行模型）。

4 实证分析

4.1 模型构建

本文以中国磷循环中畜禽养殖系统为例开展实证研究,数据来源于中国1600—2012年磷循环过程模拟研究^[4]。畜禽养殖系统（N8）的边界如图4所示,包括来自于非耕地（N2）上提供的牧草和农产品加工（N10）提供的饲料等输入量,进入农产品加工（N10）的出栏畜禽、蛋类和奶类及分别进入非耕地（N2）、内陆水体（N3）和农业种植（N7）系统里面的畜禽粪便等输出量,以及通过国际贸易（N14）输入或者输出的畜禽活体量。

图4

新窗口打开|下载原图ZIP|生成PPT
图4畜禽养殖系统边界图

Figure 4Boundary of the livestock and poultry breeding system

根据图4所示的模型系统边界,可以建立如图5所示的畜禽养殖系统磷流动态核算模型。与Spreadsheet计算方式不同的是,图5直接以图表的形式来刻画变量之间的关系,其中虚线箭头起点的变量表示为核算表达式中的自变量,终点所指的变量表示为核算表达式中的因变量,二者之间的定量关系如表1所示。图5中的变量名有两种类型,一种为椭圆,表示一个常数或计算结果;另一种为六边形,表示为一个关于时间的函数（在系统动力学里通常称为表函数或者图形函数）。前者可以表征物质流分析中确定的数值,比如含磷系数;后者可以表示随时间变化的活动水平数据,比如每年蛋类产量。在建模过程中,为了简化计算网络的复杂度,采用“影子变量”的方法实现了同一变量的多次引用（图4和图5中虚化的变量即为“影子变量”）。如果过程也需要多次引用时,可以通过构建“影子过程”来实现,比如图4中非耕地（N2）。

图5

新窗口打开|下载原图ZIP|生成PPT
图5畜禽养殖系统磷流动态核算模型图

Figure 5Phosphorus flow dynamic accounting model for the livestock and poultry breeding system

Table 1
表1
表1模型各变量计算公式
Table 1The formula of each variable in the model

变量名	计算公式
$[出栏畜禽含磷量]$	$[每年出栏畜禽含磷量] - [每年净出口畜禽含磷量]$
$[奶类含磷量]$	$[每年奶类产量] \times [奶类产品含磷系数]$
$[存栏畜禽磷的年际变化]$	$\sum^{} [畜禽 i 存栏年际变化量] \times [畜禽 i 活体含磷系数]$
$[存栏食草畜禽磷的年际变化]$	$\sum^{} [食草畜禽 i 存栏年际变化量] \times [食草畜禽 i 活体含磷系数]$
$[排入农田畜禽粪便含磷量$ ]	$\sum^{} [每年畜禽 i 粪便含磷量] \times [畜禽 i 粪便还田率]$
$[排入水体畜禽粪便含磷量]$	$\sum^{} [每年畜禽 i 粪便含磷量] \times [畜禽 i 粪便排入水体率]$
$[排入非耕地畜禽粪便含磷量]$	$\sum^{} [每年畜禽 i 粪便含磷量] \times [畜禽 i 粪便排入非耕地率]$
$[每年畜禽 i 粪便含磷量]$	$[畜禽 i 日均饲养量] \times [畜禽 i 日均排磷系数] \times 365 / 1000000$
$[每年净出口畜禽含磷量]$	$\sum^{} [畜禽 i 净出口活体数量] \times [畜禽 i 活体含磷系数]$
$[每年出栏畜禽含磷量]$	$\sum^{} [畜禽 i 出栏量] \times [畜禽 i 活体含磷系数]$
$[活体贸易含磷量]$	$[每年净出口畜禽含磷量]$
$[牧草含磷量]$	$[食草畜禽活体含磷量] + [存栏食草畜禽磷的年际变化] + [食草畜禽粪便含磷量] + [奶类含磷量]$
$[畜禽粪便含磷量]$	$\sum^{} [每年畜禽 i 粪便含磷量]$
$[粪便堆存含磷量]$	$[排入非耕地畜禽粪便含磷量]$
$[粪便排放含磷量]$	$[排入水体畜禽粪便含磷量]$
$[粪便还田含磷量]$	$[排入农田畜禽粪便含磷量]$
$[蛋类含磷量]$	$[每年蛋类产量] \times [蛋类产品含磷系数]$
$[食草畜禽活体含磷量]$	$\sum^{} [食草畜禽 i 出栏量] \times [食草畜禽 i 活体含磷系数]$
$[食草畜禽粪便含磷量]$	$\sum^{} [每年食草畜禽 i 粪便含磷量]$
$[饲料含磷量]$	$[存栏畜禽磷的年际变化] + [每年出栏畜禽含磷量] + [畜禽粪便含磷量] + [蛋类含磷量] - [牧草含磷量]$

注：本文中畜禽包括猪、牛、羊、马、驴、骡、骆驼、家禽和兔,其中食草畜禽包括牛、羊、马、驴、骡和兔。

新窗口打开|下载CSV

图5中的磷流核算包括三方面内容：①辅助变量核算;②输出磷流核算;③输入磷流核算。其中,辅助变量核算主要解决出栏畜禽及粪便含磷量长时序核算问题,采用的都是独立核算的方法。输出磷流核算主要解决包括出栏畜禽、蛋类、奶类、粪便堆存、粪便排放、粪便还田和畜禽活体贸易含磷量长时序核算问题,采用的核算方式既包括独立核算,也包括依附型核算。输入磷流核算主要解决牧草和饲料含磷量长时序核算问题,采用的核算方法既包括依附型核算也包括平衡计算。

表1给出了畜禽养殖系统磷流动态核算的计算公式,其中第一列为需要计算的变量名称,第二列为该变量的具体计算公式,所用到的变量或者参数（方括号包围的量）为图5中虚线箭头指向该变量的所有变量或者参数。

4.2 模型求解与分析

中国畜禽养殖系统磷流动态模拟结果如图6所示。图6为1600—2012年主要的输入磷流和输出磷流的模拟结果。在时间断面上,以2012年为例,结果表明,中国畜禽养殖系统中,最主要的输入磷流包括饲料（约为194.7万t,占输入总量的58.4%）和牧草（约为138.4万t,占输入总量的41.6%）;最主要的输出磷流包括畜禽粪便还田（约为146.7万t,占输出总量的43.7%）、畜禽粪便堆存（约为82.1万t,占输出总量的24.5%）、出栏畜禽（约为73.9万t,占输出总量的22.0%）、畜禽粪便排放（约为25.4万t,占输出总量的7.6%）等。由此可见,中国畜禽养殖系统大部分的磷进入到畜禽粪便中（约75.8%）。此外,在时间趋势上,磷的输入量中,自改革开放以来,饲料的投入逐渐赶超牧草的投入,并自1990年起完全超越牧草投入,说明工业化发展改变了磷的供应格局,支撑了畜禽养殖业的大规模发展;磷的输出量中,除粪便还田量在1999年达到峰值外,其他指标随时间整体呈现上升趋势,说明磷肥的大量使用导致畜禽粪便还田量降低,而多出来的畜禽粪便堆存在非耕地中,既造成了资源的浪费也成为潜在的污染源。排入水体中磷的模拟结果表明,在1960年以后畜禽养殖系统排入水体中的磷在持续上升,其原因一方面是由于社会经济发展导致人们对于肉类的需求不断增加,另一方面人口的增长也加剧了这一趋势（图6）。模型结果与Liu等^[4]的研究结果一致,这些结果对于发现中国畜禽养殖系统磷流变化规律具有重要意义。

图6

新窗口打开|下载原图ZIP|生成PPT
图6中国畜禽养殖系统磷流动态模拟结果图

Figure 6Dynamic simulation results of phosphorus flow in China’s livestock and poultry breeding system

4.3 模型参数影响分析

中国畜禽养殖系统磷循环过程模拟模型涉及很多参数,如畜禽粪便排入水体的比例等,这些参数有些涉及磷的管控水平,同时也存在一定的不确定性。在进行富营养化控制时,经常需要回答的一个重要问题是控制什么变量的问题^[28],这个就需要进行参数的敏感性分析^[29],识别那些对结果敏感的参数。另外,在评估结果可靠性的时候,需要考虑当参数取值为一个随机分布时结果的不确定性^[30],不确定性越大则模型可靠性越差。据此,本文分析了模型参数对结果的影响（图7）。图7a比较了两种畜禽粪便排入水体率情景下模拟结果的差异,其中模拟结果1是畜禽粪便排入水体率为10%时（默认值）的结果,模拟结果2是畜禽粪便排入水体率为30%（假设值）时的结果。以2012年为例,可以看出模拟结果2几乎是模拟结果1的2倍,可见在2012年如果将畜禽粪便排入水体率变为原来的1/3时（即减少2/3）,那么排入水体畜禽粪便含磷量就会变为约1/2（即减少约1/2）,也即敏感性系数为

S = (Δ y / y) / (Δ x / x) = 3 / 4

,说明该参数对于模拟结果并不十分敏感（S<1）。图7b为考虑不同畜禽粪便排入水体率均服从均匀分布

U (0,1)

时采用Monte Carlo模拟得到的不确定性分析的结果,其中黑点表示中位数,其余面积表示不同置信水平下的不确定区域范围。可见,在该假设条件下,50%的区域范围内,模型模拟的结果还是相对可靠的。

图7

新窗口打开|下载原图ZIP|生成PPT
图7中国畜禽养殖系统磷流动态模拟结果比较

Figure 7Comparison of the dynamic simulation results of phosphorus flow in China’s livestock and poultry breeding system

5 研究前景展望

传统的物质循环过程模拟大都是基于Spreadsheet的核算思路,这种研究思路需要采集大量活动水平数据和系数。由于忽略活动水平数据在时间序列上的自相关和活动水平数据之间的互相关,最终只能揭示物质循环过程的演化格局,而无法直接解释其内在的驱动因素。如果要弄清楚系统内部变量之间的关联性并识别其中的驱动因子,则需要通过相关性分析或者回归分析等研究方法来进一步挖掘^[31,32]。此外,关于活动水平数据,现有资料只是历史数据的呈现。如果预测未来的物质循环过程,就需要先对活动水平进行预测,而系统的演化也是系统性的演化过程,一个指标变化也会引起另一个指标的变化,比如人口的变化会导致人对食物的消费量的变化,而饮食结构的变化又会导致畜禽养殖量的变化,同时也会导致饲料生产量的变化,这时如果只是单纯考虑人口变化而忽略其他变化必然会引起系统偏差。另外,在系统演化时还会存在“反馈”现象,即一个因素的变化会引起另一个因素变化,而另一个因素变化又会反作用于前一个因素,这些问题在传统的物质流分析框架下都无法处理。本文提出复杂物质循环过程模拟方法除了解决物质循环过程模拟定量表达与求解问题外,对于以上缺陷也能够在一定程度上解决,具有广泛的研究和应用前景：

（1）物质循环过程动态模拟仿真。根据前面研究可知,传统的物质流分析方法一个很重要的假设就是将时间离散化,不考虑时间尺度（一般是年尺度）内的变化,即将时间尺度内部视为稳态。这种假设导致物质循环过程被肢解为一个一个时间片段,并且通过时间片段上输入输出的差异来作为系统库存的变化。这种假设的好处在于能够降低物质流核算过程的复杂度,缺点在于并未真正地实现物质循环过程动态模拟（如非稳态条件）。而本文提出的方法刚好能够解决非稳态模拟的问题,从而更加接近真实情况。

（2）物质循环过程模拟精度提升。在时间上,采用动态物质循环过程模拟的方法可以将时间精度无限细分,而不仅仅局限单一年尺度的模拟;在空间上,本文提出的方法能够与GIS相结合,从而构建出具有空间分辨率的物质循环过程动态模拟。

（3）从格局分析转向机理解析。本文提出的分析方法从系统内部各个要素及反映各个要素特征的指标入手,通过因果关联性分析来确定要素之间的定量关系,因而能够解决传统物质流分析方法仅仅只能开展格局分析而无法回答这样的格局是如何形成以及未来如何演化的问题,因而,本文的研究思路有助于实现物质循环过程模拟从定量描述性分析向定量解释性分析转变。

6 结论

本文提出了一种复杂物质循环过程模拟模型构建及其平台实现方法,并以中国1600—2012年畜禽养殖系统磷循环过程为例开展方法验证研究,得出以下结论：

（1）通过从时间、空间及过程3个维度系统解构复杂物质循环过程,识别复杂物质循环过程最小核算单元并对其开展流分析与动态核算,再将不同的核算单元整合成全过程所有区域的复杂物质循环网络的方法,能够实现复杂物质循环过程模拟模型的构建。通过数据库、模型库、结果库、分析库4“库”分离,并采用系统动力学中常微分方程组的求解方法,能够实现复杂物质循环过程模拟模型的求解。

（2）基于本文提出的方法,在EnVirLab平台上能够实现中国1600—2012年畜禽养殖系统磷循环过程模拟模型构建与求解。通过模型求解发现,模型重现了中国畜禽养殖系统磷循环发展历程,在改革开放以前主要的磷源为牧草,而之后则逐渐转变为饲料;大部分的磷进入到畜禽粪便中,其中2012年高达75.8%,约10.0%进入水体并对水体富营养化造成影响。通过对模型参数的敏感性分析和不确定性分析,发现模型参数对模拟结果的影响并非线性,且即便假定模型参数具有较大的不确定性范围,仍然能够得到不确定性相对较小的模拟结果,说明模型模拟具有较高稳健性。

参考文献原文顺序
文献年度倒序
文中引用次数倒序
被引期刊影响因子

[1]

成升魁, 甄霖. 资源流动研究的理论框架与决策应用
[J]. 资源科学, 2007,29(3):37-44.

[本文引用: 1]

[ Cheng S

, Zhen

. Resource flow: Theoretical framework and application for decision making
[J]. Resources Science, 2007,29(3):37-44.]

[本文引用: 1]

[2]

Chen M

, Graedel T

. A half-century of global phosphorus flows, stocks, production, consumption, recycling, and environmental impacts
[J]. Global Environmental Change, 2016,36:139-152.

[本文引用: 1]

[3]

陈敏鹏, 郭宝玲, 刘昱, 等. 磷元素物质流分析研究进展
[J]. 生态学报, 2015,35(20):6891-6900.

[本文引用: 1]

[ Chen M

, Guo B

, Liu

, et al. Research on phosphorus flow analysis: Progress and perspectives
[J]. Acta Ecologica Sinica, 2015,35(20):6891-6900.]

[本文引用: 1]

[4]

Liu

, Sheng

, Jiang S

, et al. Intensification of phosphorus cycling in China since the 1600s
[J]. Proceedings of The National Academy of Sciences, 2016,113(10):2609-2614.

[本文引用: 4]

[5]

Ulrich A

, Frossard

. On the history of a reoccurring concept: Phosphorus scarcity
[J]. Science of the Total Environment, 2014,490:694-707.

[本文引用: 1]

[6]

Walan

, Davidsson

, Johansson

, et al. Phosphate rock production and depletion: Regional disaggregated modeling and global implications
[J]. Resources, Conservation and Recycling, 2014,93:178-187.

[本文引用: 1]

[7]

中华人民共和国生态环境部. 2019中国生态环境状况公报
[EB/OL]. (2020-06-02) [2021-03-05]. http://www.mee.gov.cn/hjzl/sthjzk/zghjzkgb/202006/P020200602509464172096.pdf.

URL [本文引用: 1]

[ Ministry of Ecology and Environment of the People’s Republic of China. China’s Ecology and Environment Statement, 2019
[EB/OL]. (2020-06-02) [2021-03-05]. http://www.mee.gov.cn/hjzl/sthjzk/zghjzkgb/202006/P020200602509464172096.pdf.

URL [本文引用: 1]

[8]

Rockstr?m

, Steffen

, Noonk

, et al. A safe operating space for humanity
[J]. Nature, 2009,461(7263):472-475.

URLPMID:19779433 [本文引用: 1]

[9]

Steffen

, Richardson

, Rockstr?m

, et al. Planetary boundaries: Guiding human development on a changing planet
[J]. Science, 2015, DOI: 10.1126/science.1259855.

URLPMID:33888625 [本文引用: 1]

[10]

Yuan Z

, Jiang S

, Sheng

, et al. Human perturbation of the global phosphorus cycle: Changes and consequences
[J]. Environmental Science & Technology, 2018,52(5):2438-2450.

URLPMID:29402084 [本文引用: 1]

[11]

Ma D

, Hu S

, Chen D

, et al. The temporal evolution of anthropogenic phosphorus consumption in China and its environmental implications
[J]. Journal of Industrial Ecology, 2013,17(4):566-577.

[本文引用: 1]

[12]

郝敏, 陈伟强, 马梓洁, 等. 2000-2015 年中国铜废碎料贸易及效益风险分析
[J]. 资源科学, 2020,42(8):1515-1526.

[本文引用: 1]

[ Hao

, Chen W

, Ma Z

, et al. Benefits and risks of China’s copper waste and scrap trade during 2000-2015
[J]. Resources Science, 2020,42(8):1515-1526.]

[本文引用: 1]

[13]

栾晓玉, 刘巍, 崔兆杰, 等. 基于物质流分析的中国塑料资源代谢研究
[J]. 资源科学, 2020,42(2):372-382.

[本文引用: 1]

[ Luan X

, Liu

, Cui Z

, et al. Plastic resources metabolism in China based on material flow analysis
[J]. Resources Science, 2020,42(2):372-382.]

[本文引用: 1]

[14]

吴明, 姜国强, 贾冯睿, 等. 基于物质流和生命周期分析的石油行业碳排放
[J]. 资源科学, 2018,40(6):1287-1296.

[本文引用: 1]

[ Wu

, Jiang G

, Jia F

, et al. Carbon emissions from the petroleum industry based on the analysis of material flow and life cycle
[J]. Resources Science, 2018,40(6):1287-1296.]

[本文引用: 1]

[15]

沈镭, 钟帅, 胡纾寒. 全球变化下资源利用的挑战与展望
[J]. 资源科学, 2018,40(1):1-10.

[本文引用: 1]

[ Shen

, Zhong

, Hu S

. Resource utilization under global change: Challenges and outlook
[J]. Resources Science, 2018,40(1):1-10.]

[本文引用: 1]

[16]

Brunner P

, Rechberger

. Practical Handbook of Material Flow Analysis[M]. Boca Raton: CRC Press, 2004.

[本文引用: 1]

[17]

Chowdhury R

, Moore G

, Weatherley A

, et al. A novel substance flow analysis model for analysing multi-year phosphorus flow at the regional scale
[J]. Science of the Total Environment, 2016,572:1269-1280.

[本文引用: 1]

[18]

Yuan Z

, Shi J

, Wu H

, et al. Understanding the anthropogenic phosphorus pathway with substance flow analysis at the city level
[J]. Journal of Environmental Management, 2011,92(8):2021-2028.

URLPMID:21489683 [本文引用: 1]

[19]

Yuan Z

, Wu H

, He X

, et al. A bottom-up model for quantifying anthropogenic phosphorus cycles in watersheds
[J]. Journal of Cleaner Production, 2014,84:502-508.

[本文引用: 1]

[20]

Yuan Z

, Sun

, Bi

, et al. Phosphorus flow analysis of the socioeconomic ecosystem of Shucheng County, China
[J]. Ecological Applications, 2011,21(7):2822-2832.

URLPMID:22073662 [本文引用: 1]

[21]

王倩楠. 区域与企业的碳流图分析方法研究与应用
[D]. 武汉: 华中科技大学, 2018.

[本文引用: 1]

[ Wang Q

. Research and Application of Carbon Flow Diagram Analysis Method for Region and Enterprise
[D]. Wuhan: Huazhong University of Science and Technology, 2018.]

[本文引用: 1]

[22]

Cencic

, Rechberger

. Material flow analysis with Software STAN
[J]. Journal of Environmental Engineering and Management, 2008,18:3-7.

[本文引用: 1]

[23]

Rothwell S

, Doody D

, Johnston

, et al. Phosphorus stocks and flows in an intensive livestock dominated food system
[J]. Resources, Conservation and Recycling, 2020, DOI: 10.1016/j.resconrec.2020.105065.

[本文引用: 1]

[24]

Tanzer

, Zoboli

, Zessner

, et al. Filling two needs with one deed: Potentials to simultaneously improve phosphorus and nitrogen management in Austria as an example for coupled resource management systems
[J]. Science of the Total Environment, 2018,640:894-907.

[本文引用: 1]

[25]

Cooper

, Carliell-Marquet

. A substance flow analysis of phosphorus in the UK food production and consumption system
[J]. Resources, Conservation and Recycling, 2013,74:82-100.

[本文引用: 1]

[26]

Sheng

. sfc: Substance Flow Computation
[EB/OL]. (2016-08-25) [2021-03-05]. https://cran.r-project.org/web/packages/sfc/.

URL [本文引用: 1]

[27]

Treadwell J

, Clark O

, Bennett E

. Dynamic simulation of phosphorus flows through Montreal’s food and waste systems
[J]. Resources, Conservation and Recycling, 2018,131:122-133.

[本文引用: 1]

[28]

Paerl H

, Xu

, Mccarthy M

, et al. Controlling harmful cyanobacterial blooms in a hyper-eutrophic lake (Lake Taihu, China): The need for a dual nutrient (N & P) management strategy
[J]. Water Research, 2011,45(5):1973-1983.

URLPMID:20934736 [本文引用: 1]

[29]

D?ubur

, Buchner

, Laner

. Evaluating the use of global sensitivity analysis in dynamic MFA
[J]. Journal of Industrial Ecology, 2017,21(5):1212-1225.

[本文引用: 1]

[30]

Laner

, Rechberger

, Astrup

. Systematic evaluation of uncertainty in material flow analysis
[J]. Journal of Industrial Ecology, 2014,18(6):859-870.

[本文引用: 1]

[31]

Jiang S

, Hua

, Sheng

, et al. Phosphorus footprint in China over the 1961-2050 period: Historical perspective and future prospect
[J]. Science of the Total Environment, 2019,650:687-695.

[本文引用: 1]

[32]

Van Vuuren D

, Bouwman A

, Beusen A H

. Phosphorus demand for the 1970-2100 period: A scenario analysis of resource depletion
[J]. Global Environmental Change, 2010,20(3):428-439.

[本文引用: 1]