📄 fm_sae2.m
字号:
function fm_sae2(flag)
% FM_SAE2 defines a subtransmission area equivalent
% with two loads and two LTCs. This model
% uses two state variables for the tap ratios
% of the LTCs.
%
% FM_SAE2(FLAG)
% FLAG = 0 initialization
% FLAG = 1 algebraic equations
% FLAG = 2 algebraic Jacobians
% FLAG = 3 differential equations
% FLAG = 4 state matrix
% FLAG = 5 non-windup limits
%
%see also FM_SAE1 and FM_SAE3
%
%Author: Federico Milano
%Date: 11-Nov-2002
%Version: 1.0.0
%
%E-mail: fmilano@thunderbox.uwaterloo.ca
%Web-site: http://thunderbox.uwaterloo.ca/~fmilano
%
% Copyright (C) 2002-2006 Federico Milano
global Bus DAE SAE2
for i = 1:SAE2.n
if isempty(SAE2.m1)
m1 = 1;
m2 = 1;
else
m1 = DAE.x(SAE2.m1(i));
m2 = DAE.x(SAE2.m2(i));
end
ha = Bus.int(round(SAE2.con(i,1)));
hb = Bus.int(round(SAE2.con(i,2)));
xt1 = SAE2.con(i,6);
xt2 = SAE2.con(i,7);
xa0 = SAE2.con(i,8);
xb0 = SAE2.con(i,9);
ab0 = xb0;
a1 = SAE2.con(i,10);
b1 = SAE2.con(i,11);
a2 = SAE2.con(i,12);
b2 = SAE2.con(i,13);
h1 = SAE2.con(i,14);
k1 = SAE2.con(i,15);
vrif1 = SAE2.con(i,16);
h2 = SAE2.con(i,17);
k2 = SAE2.con(i,18);
vrif2 = SAE2.con(i,19);
mmax1 = SAE2.con(i,20);
mmin1 = SAE2.con(i,21);
mmax2 = SAE2.con(i,22);
mmin2 = SAE2.con(i,23);
va = DAE.V(ha);
vb = DAE.V(hb);
delta = DAE.a(ha);
theta = DAE.a(hb);
if isempty(SAE2.m1)
hm1 = 1;
hm2 = 1;
else
hm1 = SAE2.m1(i);
hm2 = SAE2.m2(i);
end
switch flag
case 0 % initialization
if length(SAE2.con(1,:)) > 20
xeq = zeros(4,1);
A = [1 0 1 0; 1 0 0 1; 0 1 1 0; 0 1 0 1; 1 1 0 0; 0 0 1 1];
x1 = SAE2.con(i,26) + SAE2.con(i,24);
x2 = SAE2.con(i,26) + SAE2.con(i,25) + SAE2.con(i,28);
x3 = SAE2.con(i,27) + SAE2.con(i,28) + SAE2.con(i,24);
x4 = SAE2.con(i,27) + SAE2.con(i,25);
x5 = SAE2.con(i,26) + SAE2.con(i,27) + SAE2.con(i,28);
x6 = SAE2.con(i,24) + SAE2.con(i,25) + SAE2.con(i,28);
xreal = [x1; x2; x3; x4; x5; x6];
xeq = A\xreal;
SAE2.con(i,8) = xeq(1);
SAE2.con(i,9) = xeq(2);
SAE2.con(i,6) = xeq(3);
SAE2.con(i,7) = xeq(4);
end
case 1
% Pa :
t1 = va*va;
t2 = xa0*xa0;
t4 = t1/t2;
t5 = vb*vb;
t6 = xb0*xb0;
t8 = t5/t6;
t9 = va*vb;
t10 = -delta+theta;
t11 = cos(t10);
t12 = 1/xa0;
t13 = t11*t12;
t19 = 1/xt1;
t20 = m1*m1;
t24 = 1/(t19+a1/t20);
t29 = 1/xt2;
t30 = m2*m2;
t34 = 1/(t29+a2/t30);
t37 = -b1/m1*t19*t24-b2/m2*t29*t34;
t38 = t37*t37;
t40 = sqrt(t4+t8+2.0*t9*t13/ab0-t38);
t42 = 1/xb0;
t43 = xt1*xt1;
t46 = xt2*xt2;
Pa = va*t40/(t19+t29+t12+t42-1/t43*t24-1/t46*t34)*t12*(t37*(-va*t12-vb* ...
t42*t11)-t40*vb*t42*sin(t10))/(t4+t8+2.0*t9*t13*t42);
% C(Qa,optimized);
t1 = va*va;
t2 = 1/xa0;
t4 = xa0*xa0;
t6 = t1/t4;
t7 = vb*vb;
t8 = xb0*xb0;
t10 = t7/t8;
t11 = va*vb;
t12 = delta-theta;
t13 = cos(t12);
t14 = t13*t2;
t20 = 1/xt1;
t21 = m1*m1;
t25 = 1/(t20+a1/t21);
t30 = 1/xt2;
t31 = m2*m2;
t35 = 1/(t30+a2/t31);
t38 = -b1/m1*t20*t25-b2/m2*t30*t35;
t39 = t38*t38;
t41 = sqrt(t6+t10+2.0*t11*t14/ab0-t39);
t43 = 1/xb0;
t44 = xt1*xt1;
t47 = xt2*xt2;
Qa = t1*t2-va*t41/(t20+t30+t2+t43-1/t44*t25-1/t47*t35)*t2*(t38*vb*t43* ...
sin(t12)+t41*(va*t2+vb*t43*t13))/(t6+t10+2.0*t11*t14*t43);
% C(Pb,optimized);
t1 = va*va;
t2 = xa0*xa0;
t4 = t1/t2;
t5 = vb*vb;
t6 = xb0*xb0;
t8 = t5/t6;
t9 = va*vb;
t10 = -delta+theta;
t11 = cos(t10);
t12 = 1/xa0;
t13 = t11*t12;
t19 = 1/xt1;
t20 = m1*m1;
t24 = 1/(t19+a1/t20);
t29 = 1/xt2;
t30 = m2*m2;
t34 = 1/(t29+a2/t30);
t37 = -b1/m1*t19*t24-b2/m2*t29*t34;
t38 = t37*t37;
t40 = sqrt(t4+t8+2.0*t9*t13/ab0-t38);
t42 = 1/xb0;
t43 = xt1*xt1;
t46 = xt2*xt2;
Pb = vb*t40/(t19+t29+t12+t42-1/t43*t24-1/t46*t34)*t42*(t37*(-vb*t42-va* ...
t12*t11)+t40*va*t12*sin(t10))/(t4+t8+2.0*t9*t13*t42);
% C(Qb,optimized);
t1 = vb*vb;
t2 = 1/xb0;
t4 = va*va;
t5 = xa0*xa0;
t7 = t4/t5;
t8 = xb0*xb0;
t10 = t1/t8;
t11 = va*vb;
t12 = delta-theta;
t13 = cos(t12);
t14 = 1/xa0;
t15 = t13*t14;
t21 = 1/xt1;
t22 = m1*m1;
t26 = 1/(t21+a1/t22);
t31 = 1/xt2;
t32 = m2*m2;
t36 = 1/(t31+a2/t32);
t39 = -b1/m1*t21*t26-b2/m2*t31*t36;
t40 = t39*t39;
t42 = sqrt(t7+t10+2.0*t11*t15/ab0-t40);
t44 = xt1*xt1;
t47 = xt2*xt2;
Qb = t1*t2-vb*t42/(t21+t31+t14+t2-1/t44*t26-1/t47*t36)*t2*(-t39*va*t14* ...
sin(t12)+t42*(vb*t2+va*t14*t13))/(t7+t10+2.0*t11*t15*t2);
DAE.gp(ha) = Pa + DAE.gp(ha);
DAE.gq(ha) = Qa + DAE.gq(ha);
DAE.gp(hb) = Pb + DAE.gp(hb);
DAE.gq(hb) = Qb + DAE.gq(hb);
% Calcolo dei termini di Jlfv
case 2
% C(diff(Pa,va),optimized);
t1 = va*va;
t2 = xa0*xa0;
t3 = 1/t2;
t4 = t1*t3;
t5 = vb*vb;
t6 = xb0*xb0;
t8 = t5/t6;
t9 = va*vb;
t10 = -delta+theta;
t11 = cos(t10);
t12 = 1/xa0;
t13 = t11*t12;
t14 = 1/ab0;
t19 = 1/xt1;
t20 = m1*m1;
t24 = 1/(t19+a1/t20);
t29 = 1/xt2;
t30 = m2*m2;
t34 = 1/(t29+a2/t30);
t37 = -b1/m1*t19*t24-b2/m2*t29*t34;
t38 = t37*t37;
t40 = sqrt(t4+t8+2.0*t9*t13*t14-t38);
t41 = 1/xb0;
t42 = xt1*xt1;
t45 = xt2*xt2;
t49 = 1/(t19+t29+t12+t41-1/t42*t24-1/t45*t34);
t58 = t41*sin(t10);
t61 = t12*(t37*(-va*t12-vb*t41*t11)-t40*vb*t58);
t64 = t4+t8+2.0*t9*t13*t41;
t65 = 1/t64;
t68 = 1/t40;
t71 = va*t3;
t72 = vb*t11;
t75 = 2.0*t71+2.0*t72*t12*t14;
t80 = va*t40*t49;
t89 = t64*t64;
DAE.J12(ha,ha) = DAE.J12(ha,ha) + t40*t49*t61*t65+va*t68*t49*t61*t65*t75/2+t80*t12*(-t37*t12-t68*vb* ...
t58*t75/2)*t65-t80*t61/t89*(2.0*t71+2.0*t72*t12*t41);
% C(diff(Pa,vb),optimized);
t1 = va*va;
t2 = xa0*xa0;
t4 = t1/t2;
t5 = vb*vb;
t6 = xb0*xb0;
t7 = 1/t6;
t8 = t5*t7;
t9 = va*vb;
t10 = -delta+theta;
t11 = cos(t10);
t12 = 1/xa0;
t13 = t11*t12;
t14 = 1/ab0;
t19 = 1/xt1;
t20 = m1*m1;
t24 = 1/(t19+a1/t20);
t29 = 1/xt2;
t30 = m2*m2;
t34 = 1/(t29+a2/t30);
t37 = -b1/m1*t19*t24-b2/m2*t29*t34;
t38 = t37*t37;
t40 = sqrt(t4+t8+2.0*t9*t13*t14-t38);
t41 = 1/t40;
t43 = 1/xb0;
t44 = xt1*xt1;
t47 = xt2*xt2;
t51 = 1/(t19+t29+t12+t43-1/t44*t24-1/t47*t34);
t59 = sin(t10);
t60 = t43*t59;
t63 = t12*(t37*(-va*t12-vb*t43*t11)-t40*vb*t60);
t66 = t4+t8+2.0*t9*t13*t43;
t67 = 1/t66;
t68 = vb*t7;
t69 = va*t11;
t72 = 2.0*t68+2.0*t69*t12*t14;
t77 = va*t40*t51;
t89 = t66*t66;
DAE.J12(ha,hb) = DAE.J12(ha,hb) + va*t41*t51*t63*t67*t72/2+t77*t12*(-t37*t43*t11-t41*vb*t60*t72/2-t40 ...
*t43*t59)*t67-t77*t63/t89*(2.0*t68+2.0*t69*t12*t43);
% C(diff(Pa,delta),optimized);
t1 = va*va;
t2 = xa0*xa0;
t3 = 1/t2;
t4 = t1*t3;
t5 = vb*vb;
t6 = xb0*xb0;
t8 = t5/t6;
t9 = va*vb;
t10 = -delta+theta;
t11 = cos(t10);
t12 = 1/xa0;
t13 = t11*t12;
t14 = 1/ab0;
t19 = 1/xt1;
t20 = m1*m1;
t24 = 1/(t19+a1/t20);
t29 = 1/xt2;
t30 = m2*m2;
t34 = 1/(t29+a2/t30);
t37 = -b1/m1*t19*t24-b2/m2*t29*t34;
t38 = t37*t37;
t40 = sqrt(t4+t8+2.0*t9*t13*t14-t38);
t41 = 1/t40;
t43 = 1/xb0;
t44 = xt1*xt1;
t47 = xt2*xt2;
t51 = 1/(t19+t29+t12+t43-1/t44*t24-1/t47*t34);
t52 = t51*t3;
t59 = t40*vb;
t60 = sin(t10);
t61 = t43*t60;
t63 = t37*(-va*t12-vb*t43*t11)-t59*t61;
t66 = t4+t8+2.0*t9*t13*t43;
t67 = 1/t66;
t69 = vb*t60;
t79 = t60*t60;
t92 = t66*t66;
DAE.J11(ha,ha) = DAE.J11(ha,ha) + t1*t41*t52*t63*t67*t69*t14+va*t40*t51*t12*(-t37*vb*t61-t41*t5*t43* ...
t79*va*t12*t14+t59*t43*t11)*t67-2.0*t1*t40*t52*t63/t92*t69*t43;
% C(diff(Pa,theta),optimized);
t1 = va*va;
t2 = xa0*xa0;
t3 = 1/t2;
t4 = t1*t3;
t5 = vb*vb;
t6 = xb0*xb0;
t8 = t5/t6;
t9 = va*vb;
t10 = -delta+theta;
t11 = cos(t10);
t12 = 1/xa0;
t13 = t11*t12;
t14 = 1/ab0;
t19 = 1/xt1;
t20 = m1*m1;
t24 = 1/(t19+a1/t20);
t29 = 1/xt2;
t30 = m2*m2;
t34 = 1/(t29+a2/t30);
t37 = -b1/m1*t19*t24-b2/m2*t29*t34;
t38 = t37*t37;
t40 = sqrt(t4+t8+2.0*t9*t13*t14-t38);
t41 = 1/t40;
t43 = 1/xb0;
t44 = xt1*xt1;
t47 = xt2*xt2;
t51 = 1/(t19+t29+t12+t43-1/t44*t24-1/t47*t34);
t52 = t51*t3;
t59 = t40*vb;
t60 = sin(t10);
t61 = t43*t60;
t63 = t37*(-va*t12-vb*t43*t11)-t59*t61;
t66 = t4+t8+2.0*t9*t13*t43;
t67 = 1/t66;
t69 = vb*t60;
t79 = t60*t60;
t92 = t66*t66;
DAE.J11(ha,hb) = DAE.J11(ha,hb) -t1*t41*t52*t63*t67*t69*t14+va*t40*t51*t12*(t37*vb*t61+t41*t5*t43* ...
t79*va*t12*t14-t59*t43*t11)*t67+2.0*t1*t40*t52*t63/t92*t69*t43;
% C(diff(Qa,va),optimized);
t1 = 1/xa0;
t2 = va*t1;
t3 = va*va;
t4 = xa0*xa0;
t5 = 1/t4;
t6 = t3*t5;
t7 = vb*vb;
t8 = xb0*xb0;
t10 = t7/t8;
t11 = va*vb;
t12 = delta-theta;
t13 = cos(t12);
t14 = t13*t1;
t15 = 1/ab0;
t20 = 1/xt1;
t21 = m1*m1;
t25 = 1/(t20+a1/t21);
t30 = 1/xt2;
t31 = m2*m2;
t35 = 1/(t30+a2/t31);
t38 = -b1/m1*t20*t25-b2/m2*t30*t35;
t39 = t38*t38;
t41 = sqrt(t6+t10+2.0*t11*t14*t15-t39);
t42 = 1/xb0;
t43 = xt1*xt1;
t46 = xt2*xt2;
t50 = 1/(t20+t30+t1+t42-1/t43*t25-1/t46*t35);
t58 = t2+vb*t42*t13;
t61 = t1*(t38*vb*t42*sin(t12)+t41*t58);
t64 = t6+t10+2.0*t11*t14*t42;
t65 = 1/t64;
t68 = 1/t41;
t71 = va*t5;
t72 = vb*t13;
t75 = 2.0*t71+2.0*t72*t1*t15;
t80 = va*t41*t50;
t88 = t64*t64;
DAE.J22(ha,ha) = DAE.J22(ha,ha) + 2.0*t2-t41*t50*t61*t65-va*t68*t50*t61*t65*t75/2-t80*t1*(t68*t58*t75 ...
/2+t41*t1)*t65+t80*t61/t88*(2.0*t71+2.0*t72*t1*t42);
% C(diff(Qa,vb),optimized);
t1 = va*va;
t2 = xa0*xa0;
t4 = t1/t2;
t5 = vb*vb;
t6 = xb0*xb0;
t7 = 1/t6;
t8 = t5*t7;
t9 = va*vb;
t10 = delta-theta;
⌨️ 快捷键说明
复制代码
Ctrl + C
搜索代码
Ctrl + F
全屏模式
F11
切换主题
Ctrl + Shift + D
显示快捷键
?
增大字号
Ctrl + =
减小字号
Ctrl + -