《经济控制论》课程研究报告
姓 名 学 号 完成日期 成 绩 联系方式
一、引言:
本次课程研究旨在学习经济系统的基本控制设计与优化方法。
应用程序语言判断系统的可控性,可观性,稳定性,研究系统的极点配置及稳定区域的求解。
二、数据资料:
1、宏观经济各个变量总量平衡统计数据表
年度
国民生产总值 (亿元) Y (k )
全社会固定资产投资
(亿元) I(k)
居民消费支出
亿元) C(k)
财政支出 (亿元) G (k )
总量平衡 (I+C+G)/Y
1986 1987 1988 1989 1990 1991 1992 1993 1994 1995 1996 10201.4 11954.9 14922.3 16917.8 18598.4 21662.5 26651.9 34560.5 46670.0 57494.9 66850.5 3120.6 3791.7 4753.8 4410.4 4517.0 5594.5 8080.1 13072.3 17042.1 20019.3 22913.5 5175.0 5961.2 7633.1 8523.5 9113.2 10315.9 12459.3 15682.4 20809.8 26944.5 32152.3 2204.91 2362.18 2491.21 2823.78 3083.59 3386.62 3742.20 4642.30 5792.62 6823.72 7937.55 1.0293 1.0134 0.9970 0.9314 0.8987 0.8908 0.9111 0.9663 0.9352 0.9355 0.9425
数据来源:根据中国统计年鉴整理。
更多数据请登录国家统计局网站查找: http://www.stats.gov.cn/ 2、 根据表中统计数据得到如下回归方程: 消费函数 C(k) = 278.85 +0.5654Y(k-1) 投资函数 I(k) = 1338+3.3664[C(k)-C(k-1)]
政府购买 G(k) = 68.617k2 - 271.51k + 2590,k 为年份,设1986年为0,依次递增。 注:扩充数据的同学需要自行建立上述方程。
三、以政府购买为控制变量,建立动态乘数-加速数经济模型的状态方程。 由题意的 取 a=2;b=0.75;
得:
x1(k+1)=b*x1(k)+b*x2(k)+b*U(k)
x2(k+1)=(ab-a)*x1(k)+ab*x2(k)+ab*U(k) y(k)=x1(k)+x2(k)+U(k)
U(k) = 68.617k2 - 271.51k + 2590
四、判定该经济系统的可控性、客观性、稳定性。
syms a b
A=[b b;a*b-a a*b]
B=[b ;a*b] C=[1 1] D=1
function character(A,B,C,D) P2=[B A*B] det(P2) rank(P2)
if rank(P2)==2 disp('******* system controllable *******') else disp('******* system not controllable ***********') end
%判系统输出的能控性 Q2=[C*B C*A*B D] rank(Q2)
if rank(Q2)==1 disp('******** system output controllable *********') else disp('******* system output not controllable *********') end
%判定系统的能观性 R2=[C;C*A] rank(R2)
if rank(R2)==2 disp('****** system observable *******') else disp('****** system not observable ***********') end
%判定系统的稳定性: function coef=polycoef(S,s,n) S=subs(S,s,'s'); for i=1:n
constant=subs(S,'s',0); coef(n-i+1)=constant; S=S-constant; S=S/'s';
S=simple(S); end
coef=[1 coef];
function result=stability(A,B) syms lamda
char_Eq=poly(A,lamda);
char_Eq_coef=polycoef(char_Eq,lamda,2); cond1=1-char_Eq_coef(3);
cond2=1-char_Eq_coef(2)+char_Eq_coef(3); cond3=1+char_Eq_coef(2)+char_Eq_coef(3); disp('----stability condition----') disp([char(cond1) '>0 ']) disp([char(cond2) '>0' ])
disp([char(cond2) '0 ']) disp('---system balance---') I=sym(eye(2)); Xe=inv(I-A)*B; Xe=simple(Xe); disp(Xe);
result=[cond1;cond2;cond3];
function closeloop(A,B,C,D,K)
disp('****** characteristic of system *******'); character(A,B,C,D);
disp('*******openloop stability analysis ******:') openloop_cond=stability(A,B); A_cl =A-B*K;
disp('*******closeloop stability analysis ******:') closeloop_cond=stability(A_cl,B);
disp('******* given specific data a=2 ;b=0.75 ; ,judge openloop atability ******') a=2;b=0.75;
openloop_cond=subs(openloop_cond); closeloop_cond=subs(closeloop_cond); disp(openloop_cond);
judge=find((openloop_cond >0) ==zeros(3,1)); flag =length(judge); if flag
disp('----openloop is not stable ----'); else
disp('---openloop is stable ----'); end
disp(' ******** closeloop stability condition ********'); disp(closeloop_cond); 运行结果: syms a b k0 k1
>> A=[b b;a*b-a a*b];B=[b;a*b];K=[k0,k1]; >> C=[1 1] C =
1 1
>> D=1 D = 1
>> closeloop(A,B,C,D,K)
****** characteristic of system ******* P2 =
[ b, b^2+a*b^2] [ a*b, (a*b-a)*b+a^2*b^2] ans = -a*b^2 ans = 2
******* system controllable ******* Q2 =
[ b+a*b, (b+a*b-a)*b+(b+a*b)*a*b, ans = 1
******** system output controllable ********* R2 =
[ 1, 1] [ b+a*b-a, b+a*b] ans = 2
****** system observable *******
*******openloop stability analysis ******: ----stability condition---- 1-a*b>0 1+2*a*b+b>0 1+2*a*b+b0
---system balance--- b/(1-b) 0
*******closeloop stability analysis ******: ----stability condition---- 1-a*b+a*b*k1>0
1-b*k0+2*a*b-2*a*b*k1+b>0 1-b*k0+2*a*b-2*a*b*k1+b0 ---system balance--- b/(1+b*k0-b) 0
1]
******* given specific data a=2 ;b=0.75 ; ,judge openloop atability ****** -0.5000 4.7500 0.2500
----openloop is not stable ----
******** closeloop stability condition ******** -1/2+3/2*k1 19/4-3/4*k0-3*k1 1/4+3/4*k0
>> maple('solve(-1/2+3/2*k1>0,{k1})') ans =
{1/3
>> maple('solve(19/4-3/4*k0-3*k1>0,{k1})')
ans =
{k1
>> maple('solve(1/4+3/4*k0>0,{k0})') ans =
{-1/3
>> lim=[-1,6,-1,3] lim =
-1 6 -1 3 >> hold on;
>> fplot('[1/3,19/12-1/4*x]',lim) plot([-1/3 -1/3],[-1 6]) plot([-1/3 -1/3],[-1 6],'r')
未加反馈校正前: >> a=2;b=0.75;
>> A=subs(A);B=subs(B); >> x1=zeros(1,120); >> x2=zeros(1,120); >> x1(1)=278.85; >> x2(1)=1338; >> X=[x1;x2]; for i=1:120
U(i)=68.617*i*i-271.51*i+2590; end
>> for i=1:120
X(:,i+1)=A*X(:,i)+B*U(i); end
>> for i=1:120
Y(i)=C*X(:,i)+D*U(i); end
>> plot(Y);
可见闭环系统不稳定。
结构图:
六、设计使国民收入平稳增长的调控政策,并说明政策含义
分析:
由上图可见,当反馈参数的取值在稳定区域靠近左侧的位置时,系统的增长趋势比较稳定,而
当取值在中间和靠近右侧时,曲线的不稳定性增大,且有较大的扰动,考虑到这些因素,选择反
馈参数K0=4,k1=7/12.
七、结语
可以清晰的从图中得出,国民生产总值在不断的增长,经济稳步发展。 在对投资行为加以负反馈控制的前提下,采用促进消费或限制消费的经济政策均可起到改善经济系统稳定性的目的
《经济控制论》课程研究报告
姓 名 学 号 完成日期 成 绩 联系方式
一、引言:
本次课程研究旨在学习经济系统的基本控制设计与优化方法。
应用程序语言判断系统的可控性,可观性,稳定性,研究系统的极点配置及稳定区域的求解。
二、数据资料:
1、宏观经济各个变量总量平衡统计数据表
年度
国民生产总值 (亿元) Y (k )
全社会固定资产投资
(亿元) I(k)
居民消费支出
亿元) C(k)
财政支出 (亿元) G (k )
总量平衡 (I+C+G)/Y
1986 1987 1988 1989 1990 1991 1992 1993 1994 1995 1996 10201.4 11954.9 14922.3 16917.8 18598.4 21662.5 26651.9 34560.5 46670.0 57494.9 66850.5 3120.6 3791.7 4753.8 4410.4 4517.0 5594.5 8080.1 13072.3 17042.1 20019.3 22913.5 5175.0 5961.2 7633.1 8523.5 9113.2 10315.9 12459.3 15682.4 20809.8 26944.5 32152.3 2204.91 2362.18 2491.21 2823.78 3083.59 3386.62 3742.20 4642.30 5792.62 6823.72 7937.55 1.0293 1.0134 0.9970 0.9314 0.8987 0.8908 0.9111 0.9663 0.9352 0.9355 0.9425
数据来源:根据中国统计年鉴整理。
更多数据请登录国家统计局网站查找: http://www.stats.gov.cn/ 2、 根据表中统计数据得到如下回归方程: 消费函数 C(k) = 278.85 +0.5654Y(k-1) 投资函数 I(k) = 1338+3.3664[C(k)-C(k-1)]
政府购买 G(k) = 68.617k2 - 271.51k + 2590,k 为年份,设1986年为0,依次递增。 注:扩充数据的同学需要自行建立上述方程。
三、以政府购买为控制变量,建立动态乘数-加速数经济模型的状态方程。 由题意的 取 a=2;b=0.75;
得:
x1(k+1)=b*x1(k)+b*x2(k)+b*U(k)
x2(k+1)=(ab-a)*x1(k)+ab*x2(k)+ab*U(k) y(k)=x1(k)+x2(k)+U(k)
U(k) = 68.617k2 - 271.51k + 2590
四、判定该经济系统的可控性、客观性、稳定性。
syms a b
A=[b b;a*b-a a*b]
B=[b ;a*b] C=[1 1] D=1
function character(A,B,C,D) P2=[B A*B] det(P2) rank(P2)
if rank(P2)==2 disp('******* system controllable *******') else disp('******* system not controllable ***********') end
%判系统输出的能控性 Q2=[C*B C*A*B D] rank(Q2)
if rank(Q2)==1 disp('******** system output controllable *********') else disp('******* system output not controllable *********') end
%判定系统的能观性 R2=[C;C*A] rank(R2)
if rank(R2)==2 disp('****** system observable *******') else disp('****** system not observable ***********') end
%判定系统的稳定性: function coef=polycoef(S,s,n) S=subs(S,s,'s'); for i=1:n
constant=subs(S,'s',0); coef(n-i+1)=constant; S=S-constant; S=S/'s';
S=simple(S); end
coef=[1 coef];
function result=stability(A,B) syms lamda
char_Eq=poly(A,lamda);
char_Eq_coef=polycoef(char_Eq,lamda,2); cond1=1-char_Eq_coef(3);
cond2=1-char_Eq_coef(2)+char_Eq_coef(3); cond3=1+char_Eq_coef(2)+char_Eq_coef(3); disp('----stability condition----') disp([char(cond1) '>0 ']) disp([char(cond2) '>0' ])
disp([char(cond2) '0 ']) disp('---system balance---') I=sym(eye(2)); Xe=inv(I-A)*B; Xe=simple(Xe); disp(Xe);
result=[cond1;cond2;cond3];
function closeloop(A,B,C,D,K)
disp('****** characteristic of system *******'); character(A,B,C,D);
disp('*******openloop stability analysis ******:') openloop_cond=stability(A,B); A_cl =A-B*K;
disp('*******closeloop stability analysis ******:') closeloop_cond=stability(A_cl,B);
disp('******* given specific data a=2 ;b=0.75 ; ,judge openloop atability ******') a=2;b=0.75;
openloop_cond=subs(openloop_cond); closeloop_cond=subs(closeloop_cond); disp(openloop_cond);
judge=find((openloop_cond >0) ==zeros(3,1)); flag =length(judge); if flag
disp('----openloop is not stable ----'); else
disp('---openloop is stable ----'); end
disp(' ******** closeloop stability condition ********'); disp(closeloop_cond); 运行结果: syms a b k0 k1
>> A=[b b;a*b-a a*b];B=[b;a*b];K=[k0,k1]; >> C=[1 1] C =
1 1
>> D=1 D = 1
>> closeloop(A,B,C,D,K)
****** characteristic of system ******* P2 =
[ b, b^2+a*b^2] [ a*b, (a*b-a)*b+a^2*b^2] ans = -a*b^2 ans = 2
******* system controllable ******* Q2 =
[ b+a*b, (b+a*b-a)*b+(b+a*b)*a*b, ans = 1
******** system output controllable ********* R2 =
[ 1, 1] [ b+a*b-a, b+a*b] ans = 2
****** system observable *******
*******openloop stability analysis ******: ----stability condition---- 1-a*b>0 1+2*a*b+b>0 1+2*a*b+b0
---system balance--- b/(1-b) 0
*******closeloop stability analysis ******: ----stability condition---- 1-a*b+a*b*k1>0
1-b*k0+2*a*b-2*a*b*k1+b>0 1-b*k0+2*a*b-2*a*b*k1+b0 ---system balance--- b/(1+b*k0-b) 0
1]
******* given specific data a=2 ;b=0.75 ; ,judge openloop atability ****** -0.5000 4.7500 0.2500
----openloop is not stable ----
******** closeloop stability condition ******** -1/2+3/2*k1 19/4-3/4*k0-3*k1 1/4+3/4*k0
>> maple('solve(-1/2+3/2*k1>0,{k1})') ans =
{1/3
>> maple('solve(19/4-3/4*k0-3*k1>0,{k1})')
ans =
{k1
>> maple('solve(1/4+3/4*k0>0,{k0})') ans =
{-1/3
>> lim=[-1,6,-1,3] lim =
-1 6 -1 3 >> hold on;
>> fplot('[1/3,19/12-1/4*x]',lim) plot([-1/3 -1/3],[-1 6]) plot([-1/3 -1/3],[-1 6],'r')
未加反馈校正前: >> a=2;b=0.75;
>> A=subs(A);B=subs(B); >> x1=zeros(1,120); >> x2=zeros(1,120); >> x1(1)=278.85; >> x2(1)=1338; >> X=[x1;x2]; for i=1:120
U(i)=68.617*i*i-271.51*i+2590; end
>> for i=1:120
X(:,i+1)=A*X(:,i)+B*U(i); end
>> for i=1:120
Y(i)=C*X(:,i)+D*U(i); end
>> plot(Y);
可见闭环系统不稳定。
结构图:
六、设计使国民收入平稳增长的调控政策,并说明政策含义
分析:
由上图可见,当反馈参数的取值在稳定区域靠近左侧的位置时,系统的增长趋势比较稳定,而
当取值在中间和靠近右侧时,曲线的不稳定性增大,且有较大的扰动,考虑到这些因素,选择反
馈参数K0=4,k1=7/12.
七、结语
可以清晰的从图中得出,国民生产总值在不断的增长,经济稳步发展。 在对投资行为加以负反馈控制的前提下,采用促进消费或限制消费的经济政策均可起到改善经济系统稳定性的目的