系统动力学模型教学案

1、第10章系统动力学模型系统动力学模型(System Dynamic)是社会、经济、规划、军事等 许多领域进行战略研究的重要工具，如同物理实验室、化学实验室一 样，也被称之为战略研究实验室，自从问世以来，可以说是硕果累累。1 系统动力学概述2 系统动力学的基础知识3 系统动力学模型第1节系统动力学概述1.1概念系统动力学是一门分析研究复杂反馈系统动态行为的系统科学 方法，它是系统科学的一个分支，也是一门沟通自然科学和社会科学 领域的横向学科，实质上就是分析研究复杂反馈大系统的计算仿真方 法。系统动力学模型是指以系统动力学的理论与方法为指导， 建立用 以研究复杂地理系统动态行为的计算机仿真模型体系

2、，其主要含义如 下：1 系统动力学模型的理论基础是系统动力学的理论和方法；2 系统动力学模型的研究对象是复杂反馈大系统；3 系统动力学模型的研究容是社会经济系统发展的战略与决策 问题，故称之为计算机仿真法的“战略与策略实验室”；4 系统动力学模型的研究方法是计算机仿真实验法，但要有计算机仿真语言 DYNAMIC勺支持，如：PD PLUS VENSIM等的支持；5 系统动力学模型的关键任务是建立系统动力学模型体系；6 系统动力学模型的最终目的是社会经济系统中的战略与策略 决策问题计算机仿真实验结果，即坐标图象和二维报表；系统动力学模型建立的一般步骤是：明确问题，绘制因果关系图, 绘制系统动力学模

3、型流图，建立系统动力学模型，仿真实验，检验或 修改模型或参数，战略分析与决策。地理系统也是一个复杂的动态系统， 因此，许多地理学者认为应 用系统动力学进行地理研究将有极大潜力，并积极开展了区域发展， 城市发展，环境规划等方面的推广应用工作，因此，各类地理系统动 力学模型即应运而生。1.2 发展概况系统动力学是在20世纪50年代末由美国麻省理工学院史隆管理 学院教授福雷斯特（JAY.W.FORRESTERI出来的。目前，风靡全世 界，成为社会科学重要实验手段，它已广泛应用于社会经济管理科技 和生态灯各个领域。福雷斯特教授及其助手运用系统动力学方法对全 球问题，城市发展，企业管理等领域进行了卓有成

4、效的研究，接连发 表了工业动力学，城市动力学，世界动力学，增长的极限 等著作，引起了世界各国政府和科学家的普遍关注。在我国关于系统动力学方面的研究始于 1980年，后来，陆续做 了大量的工作，主要表现如下：1 ）人才培养自从1980年以来，我国非常重视系统动力学人才的培养，主要 采用“走出去，请进来”的办法。请进来就是请国外系统动力学专家 来华讲学，走出去就是派留学生，如：首批派出去的复旦大学管理学 院的王其藩教授等，另外，还多次举办了全国性的讲习班。2 ）编译编写专著组织专家编译了工业动力学，城市动力学等。编写专著有：王其藩著系统动力学，高级系统动力学；胡 玉奎著系统动力学，王洪斌著系统动力

5、学教程，贾仁安著系 统动力学教程等。3 ）引进专业软件引进的软件有：MICRO-DYNAMDYNAMAPDYNAMO ，STELLA PD PLUS等，近几年又引进的最先进实用的 VENSIM专业软件。并自 行研制了一些专用软件。4 ）新设课程新开设了系统动力学专业课程。在几十所大学的管理系或管理学 院以及科研单位的研究生开设了系统动力学课程。5 ）组织与学术会议于19年成立了全国系统动力学委员会。组建了一些专门研究 和教学。开展了许多专项研究工作。建立了国家总体系统动力学模型， 省和地区的发展战略研究系统动力学模型， 省级能源，环境预测系统 动力学模型及科技，工业，农业林业等行业发展战略研究

6、系统动力学 模型等。1986 年8月，在上海召开的“全国系统动力学学术研讨会“上，140多名代表提交了 95篇有关系统动力学理论和应用研究方面的论 文。1987年6月，在上海召开的国际学术会议上我国代表交流了29篇论文，占会议论文数的45% 1988年7月，美国圣迭戈召开了国际 学术年会，我国有十名代表参加，交流论文十多篇。1989年7月，在西德斯图加特召开的国际学术年会上，我国学者交流论文14篇，有4人参加会议。目前，在我国系统动力学已经发展成熟，并正向深入和全面应用 延伸，形成了一支强大的研究力量，发展趋势看好，有理由相信，系 统动力学必将在我国社会，经济，科技，管理和生态等领域的研究中

7、发挥更大作用。第2节系统动力学的基础知识系统动力学模型建立的基本知识，基本原理主要有：因果关系图， 模型流图及模型的组成等。现分别介绍。2.1因果关系1 因果关系因果关系是指由原因产生某结果的相互关系。 从哲学角度讲，原 因和结果是揭示客观事物的因果联系的重要哲学概念， 它们是客观事 物普遍联系和相互作用的表现形式之一。 原因是某种事物或现象，是 造成某种结果的条件；结果是原因所造成的事物或现象，是在一定阶 段上事物发展所达到的目标状态。通常用箭头线来表示，它有正因果关系和负因果关系两种， 如图9 1。Pi69原因结果+就业机会E迁入人口数I死亡率R总人口数P正因果关系：两个变量呈同方向变化趋

8、势，如：E增加，I增加;E减少，I减少。负因果关系：两个变量呈异方向变化趋势，如：R增加，P减少;R减少，P增加。2 ）因果关系环图因果关系环图是指由两个或两个以上的因果关系连接而成的闭合回路图示。它定性描述了系统中变量之间的因果关系。它有正负因果关系环图两种，如图93，图9-4所示：P169正因果关系环图：它会引起系统部活动加强。准则：若各因果关系均为正，则该环为正因果关系环；若各因果关系为负的个数是偶数时，则该环也为正因果关系环。负因果关系环图：它会引起系统部活动减弱准则：若各因果关系均为负，则该环为负因果关系环；若因果关系为负的个数是奇数，则该环为负因果关系 环。再如：生态学人口增长因果

9、关系环图，如图95,图9-6所示:p702.2 系统动力学模型流图系统动力学模型流图简称 SD流图，是指由专用符号组成用以表 示因果关系环中各个变量之间相互关系的图示。 它能表示出更多系统 结构和系统行为的信息，是建立 SD模型必不可少的环节，对建立 SD 模型起着重要作用。其专用符号主要有八个：1 ）水平变量水平变量符号是表示水平变量的积累状态的符号，它是 SD模型 中最主要的变量。它由五部分组成，即：输入速率，输出速率，流线， 变量名称及方程代码（L）,如图 所示。2）速率变量速率变量符号是表示水平变量变化速率的变量。 它能控制水平变 量的变化速度，是可控变量。它由三部分组成，即：输入信息

10、变量， 变量名称及方程代码（R）。如图 所示。3）辅助变量辅助变量符号是辅助水平变量等的变量。如图所示4）外生变量外生变量符号如图 所示。5）表函数表函数符号如图所示。6）常数常数符号如图所示。7）流线流线符号又有物质流线，信息流线，资金流线，及订货流线四种：物质流线符号是表示系统中流动着的实体，如图所示。信息流线符号是表示联接积累与流速的信息通道，如图所示。资金流线符号是表示资金，存款及货币的流向，如图所示。订货流线符号是表示订货量与需求量的流向，如图所示。8）源与沟源符号与沟符号如图 所示。2.3系统动力学模型系统动力学模型是由六种基本方程和专门的输出语句组成。其六种方程的标志符号分别为:

11、L:水平变量方程；:速率变量方程;A:辅助变量方程；N:计算初始值方程；C:赋值予常数方程；T:赋值予表函数中Y坐标值。L方程是积累方程；R, A方程是代数运算方程；C, T, N方程是提供参数值方程，并在同一次模拟中其值保持不 变。1)L方程L方程是计算水平变量积累值的方程，其一般表示形式为：L POP鬃k = P°P j + DT ? (BR 鬃k DR jk )其中，L :水平变量方程代码，表示方程性质。DT :时间间隔，即时间增量。J :表示前一刻。.K :现在时刻。丄：未来一时刻。POPj :过去一时刻人口数。POPk :现在时刻人口数。POPL :未来一时刻人口数。BR

12、jk :过去至现在该段时刻的人口出生率。DR JK :过去至现在该时刻段的人口死亡率。积累是系统部流的堆积量，它等于过去一时刻的积累加上积累 变动量，即变动增量。积累变动量是时间间隔与输入流速和输出流速 之差的乘积。2） R方程R方程是计算单位时间流量的方程，即流速或速率。其一般表 示形式为：RBRjkBRFPOPjRDJDRF POPjRBR klBRFPOPkRDR klDRF POPk其中，BRjk :过去至现在时刻的出生率，单位（人/年）；DR jk :过去至现在时刻的死亡率，单位（人/年）；BRkl :现在至未来时刻的出生率；单位（人/年）；DR kl :现在至未来时刻的死亡率，单位

13、（人/年）；brf :出生系数，单位（人/年.人）；DRF :死亡系数，单位（人/年.人）；POPj :过去时刻人口总数；POPk :现在时刻人口总数。3 ） A方程A 方程是辅助变量方程，用于对辅助变量赋值，其一般表示 形式为：A TPOP K sum（22, popQ其中，TPOPk :表示现在人口总数。SUM （22, pop k）:求和函数，表示求算现在 22个年龄组的总和。4 ) N方程N 方程是变量初始值方程，表示对变量赋初始值，起一般表示形式为：NPOP(AGEI) IPOP(AGEI)其中，POP(AGE 1):表示各年龄组人口初始值。IPOP(AGEI):是表函数，表示存储2

14、2个年龄组的初始值5 )T方程T 方程是表函数方程，表示对相应的纵坐标Y赋值。6 )C方程C方程是常数方程，表示对常数变量赋值。第3节系统动力学模型系统动力学模型应用分析的一般步骤为:1 明确问题明确的问题是：系统的围：空间围，如省区域；时间围，女口 1961年-2050 年;时间间隔，DT=1年，等等。解决途径：计算机仿真实验。数据资料：人口总数，出生率，死亡率，自然增长率等。2 明确目标人口总数变化趋势；自然增长率控制目标；出生率控制目标；死亡率控制目标等。3 绘制系统流图1 )因果关系环图主要变量清单，即列出主要变量的清单，以利于因果关系环流图 的绘制。如：总人口数，出生率，死亡率，出生

15、系数，死亡系数。很容易 绘制出下图：2 )SD模型流图在因果关系环图的基础上可得 SD模型流图如图 所示4 SD 模型的建立根据上述介绍知识和分析步骤，可得简单的省人口 SD模型如下:* POPULA YION SD MODEL OF ANHUIPOPk POP J DT (BRjkDRjk)RBR KLBRFPOPkRDR KLDRFPOPk)NPOP60000000CBRF0.005CDRF0.003SPEC DT=1/PRINT 1） POP 2） BR 3） DRPLOT POP BR DRPLOT POP说明：1 ）人口数分22个年龄组，即：1 岁，2 4 , 5 9 , 10 14

16、 ,。，95 99 , 100 岁数以上计22个年组。251等。2）增加与迁移人口在实际应用时要包括：增加人口和迁移人口。模拟或仿真运行）进入SD状态，即SD专业软件状态，如：PD PLUS VENSIM）输入程序）运行程序（1）按初始设计参数，100年总人口的变化情况；（2） 调整人口出生系数和人口死亡系数，100年总人口的变化3情况;（3）确定各个时期的人口目标，实验取得最佳的BRF制定合理的人口生育政策。6 解释分析模拟结果1 ）对比分析 按时间段纵向分析，按区域分布横向分析。2 ）效果分析 参数的灵敏度，政策的实用性和灵活性。3 ）技术分析SD的理论与方法技术，模型建立技术与技巧，SD

17、专业软件操作技术等。7 修正再模拟1 ）模型结构适合性检验（1）量纲量纲必须正确。（2）极端条件极端条件下是否仍有意义。（3）边界条件2 ）模型行为适合性检验（1）结构灵敏度（2）参数灵敏度（3）结合灵敏度3 ）模型结构与真实系统一致性4 ）模型行为与真实系统一致性现将省人口 SD模型列于以下：* POPULATION SD MODEL OF ANHUIC NAGES=22FOR AGE1 = 1，NAGES=AGAG2 AG3 AG4 AG5 AG6 AG7 AG8 AG9 AG10 AG12, AG13 AG14 AG15 AG16 AG17 AG18 AG19 AG20 AG21 AG2

18、2FOR AGE2=2，NAGESA TPOPk sum(22, POPk)L POP k(AG1) POP j(AG1) DT (BRjk DRjk(AGI) MOUTjk(AGI)POPk(AGE2) POPj(AGE2) DT (MOUT jk(AGE2 1)MOUT j k (AGE2) DR j k (AGE2)n pop(agei =ipop(ageiT IPOP(*)=867601,3881562,5772587,7431281,6740139,24887984108669，3507741，2741800，2486142，2406997, 2031682,/p>

19、043 ， 945653, 513026, 352871， 166180, 42022， 7165 881,121,29!29不详R BRkl SCLPRD(BRF( ),1,7,POPk( ),5) !积函数，即 16 20 至 49 岁计7组。T BRF=0.0457795 , 0.07161650, 0.11789292, 0.04525700.01577322, 0.00636126 , 0.00148608R DRkl(AGEI) DRF(AGEI) POPk(AGEI)T DRF=0.02216, 0.004725, 0.001 , 0.00066 , 0.000960.00156

20、,0.00152 , 0.00173,0.00213,0.00307,0.004470.00741,0.01154,0.01891,0.02899,0.04986,0.07005,0.0801,0.1751,0.8851,0.9660,0.9752! 22 个年龄组的死亡系数。R MOUT kl(AGE1) POPk(AGE1)/MT(AGE1) !表示每组每年转出人 数。T MT=1, 4, 5, 5, 5, 5, 5, 5, 5, 5, 5, 5, 5, 5, 5, 5, 5, 55, 5, 5, 1E+30! 1E+30表示最后一组转出为零。SPEC DT=1/LENGTH=100/PR

21、TPER=2/PLTPER=1/SAVPER=2PRINT 1)TPOP 2)POP 3)BR 4)DRplot tpop br drplot pop第4节系统动力学模型的实习指导4.1 Vensim的基本操作系统动力学模型是由SD专业软件来支撑仿真运行的，例如 SD专业软件 Micro-Dynamo I、H, PD PLUS Vensim等。但只有 Vensim在Windows下运行，而且是非常方便、实用的新版本SD专业软件故仅介绍Vensim的基本操作。4.1.1安装与启动(1) 安装将装有安装文件和Vensim系统文件一同COPY到硬盘，然后单击安装文件依据提示即可安装。(2) 启动安装

22、成功后，快捷方式 Vensim图标显示在桌面上，双击该图标即进入Vensim窗口。如图8-2所示。M1MLtMJWLAILJ；RIHRATI*lUbbii7号“一兰乂理、更fP£? 创ild giwtrol Digut JJJp 哼nhfeit birthsRABBI1 UFl o7 j懸7«pfcl*?nMA1HJ LOOKURj图 8-2 Vensim 主窗口4.1.2 窗口组成(1)主菜单File文件管理New Modle 建立新模型Ope n Modle打开模型Save保存模型Pri nt打印模型Check检查模型View 观察：对模型的流图进行总体观察Datase

23、ts数据集分析：Simulate运行模型Reality Check真实性检验Time Bou nds时间控制设置Graph自定义作图：用于由用户选择作图。Build模型建立：实现由其它窗口转换到模型建立窗口Co ntrol控制窗口：Variable Selecti on分析变量选择Time Axis时间轴设置Scalli ng纵坐标调整Datasets数据集选择Custom Graph自定义作图Output输出窗口：实现从其它跳转到输出窗口(2)主窗口主窗口由三个子窗口组成。Buildi ng Win dow 模型建立窗口：用于建立模型Con trol Win dow控制窗口：用于调整和控制输

24、出结果，使输出曲线等更加完美Output Window输出窗口：用于结构分析和数据集分析选择三个窗口可通过工具条：Building /Output/Control实现切换。4.1.3模型建立模型建立过程是先输入流图，然后输入模型的方程和参数。(1) 打开模型建立窗口单击“ File ”菜单，弹出子菜单；单击“ New Modle” 或“ Open Modle” 即可。(2) 输入流图 绘制流图根据分析的因果图绘制流图如图8-2所示。 输入流图采用系统的图示工具条，位于主菜单下的工具行，很容易输入流 图。主要工具及使用方法，其功能从左向右依次为：Poiny指针：用于选中、移动和改变变量和图符。选

25、择方式有三种：一是按下左键拖动一个虚线方框，松开后则方框所有变量和图符 被选中，选中后可用鼠标拖动之；二是直接单击变量或图符选中，然后可按下鼠标左键拖动，另外 若选中图符操作柄（一个小圆圈），可对图符进行放大缩小操作；三是双击变量将其选为工作台变量。Variable 变量：用于定义非状态变量，如辅助变量和常量。BoxVariable方框变量：用于创建状态变量。注意，Vensim 一般约定，状态变量第一个字母大写，常量和表 格量字母全大写，其余变量字母全小写。但用户也可以有自己的约定。Rates速率变量：用于创建速率变量。它有四部分构成，两个箭头，一个开关，两 朵表示源和漏的云和变量本身。 可通

26、过选择移动开关，移动和改变图 符形状。速率变量一般至少有一端指向一个状态变量，在创建时可是 其一端直接指到状态变量上。Arrow箭头：用于创建表示因果关系的前头，可以是直的或弯曲的。注意：创建时，先选中箭头，在因变量上单击左键，然后松开， 移动鼠标到空白处单击左键，松开后移动到果变量再单击左键。Comments 释:用于为流图增加注释，使容易看懂。注释方法相当丰富，可以是 文字，图符等。Shadow Variable再声明变量：若一个变量已在一个子块中定义，在另一个子块中又需应用时， 只要用其再声明一次即可。 编辑流图图示编辑工具条位于主窗口最下一行，计有16种工具，专用于对流图进行编辑，非常

27、方便，光标移至某个工具上时会显示提示信息。(3) 输入方程和参数输入方程的具体过程是：打开方程编辑对话框打开该对话框的方法有两种：一是鼠标右击流图某一变量，弹出变量设置对话框，再单击“ Equation ”项即可打开方程编辑对话框如图 8-3所示。n J*亡*1»:FUbbii PopuhUonFrlc . E4H View 0»ta£ct3 Grattht Build Control Oulnul He tn图8-3方程编辑对话框二是在建立模型窗口选择 Equation Editor，然后单击某一变量,即进入该对话框输入方程和参数采用该对话框很容易实现方程输入，

28、主要功能是：方程编辑：可选择变量、函数、数字和运算符构成方程。In put输入变量：若流图正确，因变量均包含于此列表框中，并构成输入变量，在 方程编辑时只要从此选择即可。Function 函数：Vensim提供的所有函数均在此列表框中，只要选择即可进入方 程编辑。Lookup Table 查找表:若方程左边包含表函数即可选此项。Un its单位：方程右边必须有正确单位。若是一个无量纲，可不加入或输入dmnlComments 释：可对方程注释，以增加可读性。Erros错误信息：在方程书写时，会随时提示书写是否正确。4.1.4运行模型模型的方程和参数输入完毕，即可运行模型。具体过程是：首先进入Ve

29、nsim主窗口，然后单击工具条上的“运行”钮，即进行仿真运行。4.1.5输出运行结果输出运行结果是采用数据分析的工具条实现的。若单击“变量曲线”钮，即显示变量的坐标曲线图；若单击“数据表”钮，即显示仿真结果数据表。4.2系统动力学模型的实习指导421 实习目的1、巩固系统动力学仿真实验的基本原理和方法步骤；2、掌握系统动力学模型的建立、输入、仿真、输出及应用方法 与技巧；3、求取系统动力学模型的仿真结果并应用于科学决策；4、掌握系统动力学模型的变换应用方法。4.2.2 实习容实习容为我国水资源可持续发展系统动力学模型。2、主要变量说明TWQ水资源可持续发展供水量MUTWR可用总水量RUWR重复

30、利用水率SDC可持续发展开采水资源因子LWC生活配水系数PQ人口数PIQ人口出生率IWRQ工业用水需求WTPC水技术人员因子WTIC水技术进步因子PWDQ污水排放量PWLR污水比率PWAT污水净化时间TWPWC供水导致污水因子PWLS污水率参考值TWC第三产业配水系数 TPFC农业供水事业投资比率AWC农业配水系数TTND第三产业供水紧程度IWC工业配水系数IPAIV工业产值实际增长速度ITW工业供水量ITND工业供水紧程度IP工业产值LRWQ生活需水量ITPFC工业供水事业投资比率WTP水技术人员数PDQ人口死亡率ATW农业供水量ATND农业供水紧程度APAIV农业产值实际增长速度AP农业产

31、值TP第三产业产值TTWI第三产业供水事业投资TTPFC第三产业供水事业投资比率TWRQ第三产业用水需求AWRQ农业用水需求TTW第三产业供水量TPAIV第三产业产值实际增长速度LTND生活供水紧程度PDR人口死亡率控制因子IWIR工业用水增长率3、模型该模型的主要方程是:TWQ二MUTWQ*TWCC*(1+RUWR+SDC+WTIC*WTP)*ATWC*TTWC*ITWCPDQ二PDR*PQ*PWDQ/LTNDITWI=IFCI*ITPFCIP=INTEG(IPR, 5.251e+011)IWRIV=IP*IWAPCIPR=IP*IPAIV*IPAIVCIFCI=IFPAC*IPIWRQ=I

32、NTEG(IWIR, 3.2E+)PWLGR=TWPWC*TWQPWLAR二PWDQ/PWATPWDQ=INTEG (+PWLGR-PWLAR, 3.487E+)PQ=INTEG(+PIQ-PDQ, 1.093E+009)LWRQ=INTEG(LWIR, 1.97E+)LWRIV=LWAPC*PQLWIR=LWRIV*LWRQ(注：三个产业的方程中以工业为例)4.2.3实习过程以我国水资源可持续发展系统动力学模型为实习容，其具体实习过程是：1、绘制流图在系统分析的基础上，先绘制该模型的因果图，然后再绘制出流图。2、编写模型的方程根据流图编写出模型所有方程及参数值的确定。3、输入流图采用图示工具

33、，将绘制的流图输入计算机。4、输入方程和参数采用模型建立窗口提供的功能，将模型的方程和辅助方程输入计 算机。5、运行模型模型方程及辅助方程输入完毕，即可执行运行命令进行模型模 拟。6、输出运行结果单击数据集分析工具条中的“变量曲线”钮，即显示变量曲线图5一凸-s/fu看趋聿$2石歩.1 覺载1=SXO 壬*華7、仿真模拟运行仿真模拟运行就是修改参数进行多方案试验模拟运行，以选择最佳方案，供科学决策参考。424 SD模型的变换应用再次通过人口系统动力学模型的建立、输入、运行、模拟等整个 过程，使进一步掌握SD模型的应用方法与技巧。1、准备工作(1) 绘制流图人口模型流图如图8-5图8-5 人口

34、SD模型流图(2) 编写人口 SD模型人口 SD的方程及辅助方程如下:*POPULATIONMODELLP.K二P.J+DT*(BR .J K-DR .JK)NP=20130000NOTEP- POPULATION BR BIRTH RATE DR DEATHRATERBR.KL=BW.K*NBR*CC.KARBR.K=NBR*CC.KNOTERBRREAL BIRTHRATEANAL.K=RBR.K*CBPCCBP=30NOTENAL- NATALITY FOR LIFENOTECBP- CHILD-BEARING PERIODNOTEBV BIRTH WOMENNOTECC CONTROL

35、 COEFFICIENTNOTENBR- NORMAL BIRTH RATE OF BIRTH WOMENCNBR=0.25NTIME=1964ABW.K=P.K*PRPCPRP=0.2NOTEPRP- THE PROPORTION OF BW IN PRDR.KL=P.K*AMTCAMT=0.0061NOTEAMT- AVERAGE MORTALIT YAPD.K二P.KZ AREACAREA=473250NOTEPD- POPULATION DENSITYACC.K二TABHL(TCC, PD.K, 40,80,10)TTCC=T .85 /.6 /.2 /.2NOTETCC- CONTROL COEFFICIENT TABLENOTECONTROL STATEMENTSPECDT=1/ LENGTH=2004PRTPER=2PLTPER=2PRINTP,BR,DR,PD,BW,CC,RBR,NALPLOTP=P/ BW二B(4000000,10000000)/ BR=R, DR=D4PD=8*RUNRUN1(3) 运行结果参考运