-
Notifications
You must be signed in to change notification settings - Fork 4
/
Copy pathhFunctionVectorized.m
72 lines (56 loc) · 2.25 KB
/
hFunctionVectorized.m
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
function [h1,h2,h3,h4,h5,h6] = hFunctionVectorized(delta,e,v,theta,...
pg,qg)
%HFUNCTIONVECTORIZED computes h(x,a,u); i.e., algebraic equations in DAE.
% [h1,h2,h3,h4,h5,h6] = hFunction(delta,e,theta,V,...
% pg,qg,xq_vec,xprime_vec) computes the algebraic equations of
% function h based on equations (2a)--(3b)
%
% Description of outputs:
% 1. h1: size(G,1) right-hand side of equation (2a)
% 2. h2: size(G,1) right-hand side of equation (2b)
% 3. h3: size(G,1) right-hand side of equation (2c)
% 4. h4: size(G,1) right-hand side of equation (2d)
% 5. h5: size(L,1) right-hand side of equation (3a)
% 6. h6: size(L,1) right-hand side of equation (3b)
%
% Description of inputs:
% 1. delta: size(G,1) instantaneous internal angle of generator in radians
% 2. e: size(G,1) instantaneous electromotive force
% 3. theta: size(N,1) instantaneous terminal angle of nodes in radians
% 4. V: size(N,1) instantaneous terminal voltage magnitude of nodes in
% pu volts
% 5. pg: size(G,1) instantaneous real power generated by the generator
% 6. qg: size(G,1) instantaneous reactive power generated by the generator
% 7. xq_vec: size(G,1) vector of quadrature axis synchronous reactance (pu)
% 8. xprime_vec: size(G,1) vector of direct axis transient reactance (pu)
%
% See also hFunction
% h1--h4 size(G,1)
% h5,h6 size(L,1)
global G L gen_set load_set Gmat Bmat Cg...
global xq_vec xprime_vec
% constants:
Xprime=diag(xprime_vec);
Xq=diag(xq_vec);
% vector variables:
thetag=theta(gen_set);
Vg=v(gen_set);
% matrix variables
E=diag(e);
V_g=diag(Vg);
Vmat=diag(v);
cosMat=diag( cos(theta));
sinMat=diag( sin(theta));
h1=-pg+ inv(Xprime)*E*V_g*sin(delta-thetag)+...
(1/2)*inv(Xq)*inv(Xprime)*(Xprime-Xq)*V_g*V_g*sin(2*(delta-thetag));
h2=-qg+inv(Xprime)*E*V_g*cos(delta-thetag)-...
(1/2)*inv(Xq)*inv(Xprime)*(Xprime+Xq)*V_g*Vg+ ...
(1/2)*inv(Xq)*inv(Xprime)*(Xprime-Xq)*V_g*V_g*cos(2*(delta-thetag));
h35=Vmat*cosMat*Gmat*Vmat*cos(theta)-Vmat*cosMat*Bmat*Vmat*sin(theta)...
+Vmat*sinMat*Bmat*Vmat*cos(theta)+Vmat*sinMat*Gmat*Vmat*sin(theta)-Cg*pg;
h3=h35(gen_set);
h5=h35(load_set);
h46=Vmat*sinMat*Gmat*Vmat*cos(theta)-Vmat*sinMat*Bmat*Vmat*sin(theta)...
-Vmat*cosMat*Bmat*Vmat*cos(theta)- Vmat*cosMat*Gmat*Vmat*sin(theta)-Cg*qg;
h4=h46(gen_set);
h6=h46(load_set);