How can the calculation time of the system matrix be reduced in this
- To: mathgroup at smc.vnet.net
- Subject: [mg84619] How can the calculation time of the system matrix be reduced in this
- From: lightnation <lightnation at naver.com>
- Date: Mon, 7 Jan 2008 02:38:38 -0500 (EST)
The tagrget is to get the system matrix in the form of:
x' = A.x
The problem is,
A is expressed in the alegraic form, not in the form of numerical
form.
I want to reduce the calculation time of this system matrix.
In this file, Asyspssonly is the system matrix A.
The main obstacle for the calculation time is caused by Inverse[D75]
in the section of 'Matrix manipulation 3'
=================================================
Equations for powerflow and dynamic stability
Line Data
In[114]:= Off[General::"spelL"]
Off[General::"spell1"]
linedata = {{1, 4, 1, I 0.0576, 0, 0}, {2, 7, 1, I 0.0625, 0, 0}, {3,
9, 1,
I 0.0586, 0, 0}, {4, 5, 1, 0.01 + I 0.085, I 0.088, I 0.088}, {4,
6, 1,
0.017 + I 0.092, I 0.079, I 0.079}, {5, 7, 1, 0.032 + I 0.161, I
0.153,
I 0.153}, {6, 9, 1, 0.039 + I 0.17, 0.179 I, 0.179 I}, {7, 8, 1,
0.0085 + I 0.072, I 0.0745, I 0.0745}, {8, 9, 1, 0.0119 + I
0.1008,
I 0.1045, I 0.1045}};
Y matrix with line data
In[16]:= FY1[linedata_] :=
Module[{},
index = Max[
Flatten[Table[{linedata[[i, 1]], linedata[[i, 2]]}, {i, 1,
Length[linedata]}]]];
trn = Max[Union[Table[linedata[[i, 3]], {i, 1,
Length[linedata]}]]];
Do[Subscript[linedata, i] = {}, {i, 1, trn}];
For[i = 1, i <= Length[linedata],
For[q = 1, q <= trn,
If[linedata[[i, 3]] == q,
AppendTo[Subscript[linedata, q], linedata[[i]]],]; q++]; i++];
Ymat1 = ZeroMatrix[index, index];
Ymat2 = ZeroMatrix[index, index];
FY2 = ZeroMatrix[index, index];
For[q = 1, q <= trn,
Map[(Ymat1[[#[[1]], #[[2]]]] = #[[4]]^-1) &, Subscript[linedata,
q]];
Map[(Ymat1[[#[[2]], #[[1]]]] = #[[4]]^-1) &, Subscript[linedata,
q]];
Map[(Ymat2[[#[[1]], #[[2]]]] = #[[5]]) &, Subscript[linedata, q]];
Map[(Ymat2[[#[[2]], #[[1]]]] = #[[6]]) &, Subscript[linedata, q]];
For[i = 1, i <= Length[Ymat2],
For[j = 1, j <= Length[Ymat2],
If[i == j,
FY2[[i, j]] =
FY2[[i, j]] +
Sum[Ymat1[[i, k]] + Ymat2[[i, k]], {k, 1, Length[Ymat1]}],
FY2[[i, j]] = FY2[[i, j]] - Ymat1[[i, j]]]; j++]; i++];
Ymat1 = ZeroMatrix[index, index];
Ymat2 = ZeroMatrix[index, index]; q++]; Return[FY2]];
Power Flow Function
In[17]:= PowerFlow[BUSdat1_, linedata_, SV_] :=
Module[{mismatch1, mismatchcal}, BUSN = Length[BUSdat1];
GBus = {};
DBus = {};
FY = FY1[linedata];
For[i = 1, i <= BUSN,
If[BUSdat1[[i, 2]] == 2, AppendTo[GBus, BUSdat1[[i]]]]; i++];
For[i = 1, i <= BUSN,
If[BUSdat1[[i, 2]] == 1, AppendTo[DBus, BUSdat1[[i]]]]; i++];
GenIndex = Table[GBus[[i, 1]], {i, 1, Length[GBus]}];
LoadIndex = Table[DBus[[i, 1]], {i, 1, Length[DBus]}];
Do[Subscript[Peq, i] = Subscript[V, i] \!\(
\*UnderoverscriptBox[\(\[Sum]\), \(k = 1\), \(BUSN\)]
SubscriptBox[\(V\), \(k\)] \((Re[FY[\([\)\(i, k\)\(]\)]] Cos[
\*SubscriptBox[\(\[Theta]\), \(i\)] -
\*SubscriptBox[\(\[Theta]\), \(k\)]] + Im[FY[\([\)\(i, k\)\(]\)]]\
Sin[
\*SubscriptBox[\(\[Theta]\), \(i\)] -
\*SubscriptBox[\(\[Theta]\), \(k\)]])\)\);
Subscript[Qeq, i] = Subscript[V, i] \!\(
\*UnderoverscriptBox[\(\[Sum]\), \(k = 1\), \(BUSN\)]
SubscriptBox[\(V\), \(k\)] \((Re[FY[\([\)\(i, k\)\(]\)]] Sin[
\*SubscriptBox[\(\[Theta]\), \(i\)] -
\*SubscriptBox[\(\[Theta]\), \(k\)]] - Im[FY[\([\)\(i, k\)\(]\)]]\
Cos[
\*SubscriptBox[\(\[Theta]\), \(i\)] -
\*SubscriptBox[\(\[Theta]\), \(k\)]])\)\), {i, 1, BUSN}];
SVA = Thread[{Subscript[\[Theta], SV], Subscript[V,
SV]} -> {BUSdat1[[SV, 6]], BUSdat1[[SV, 5]]}];
var = Flatten[{Map[Subscript[\[Theta], #] &, Complement[GenIndex,
{SV}]],
Map[Subscript[\[Theta], #] &, LoadIndex],
Map[Subscript[V, #] &, LoadIndex]}];
KV = Thread[
Map[(Subscript[V, #]) &, Complement[GenIndex, {SV}]] ->
Map[ (BUSdat1[[#, 5]]) &, Complement[GenIndex, {SV}]]];
varcal0 =
Flatten[{Map[BUSdat1[[#, 6]] &, Complement[GenIndex, {SV}]],
Map[(BUSdat1[[#, 6]]) &, LoadIndex],
Map[(BUSdat1[[#, 5]]) &, LoadIndex]}];
PQS0 = Flatten[{Map[BUSdat1[[#, 3]] &, Complement[GenIndex,
{SV}]],
Map[BUSdat1[[#, 3]] &, LoadIndex], Map[BUSdat1[[#, 4]] &,
LoadIndex]}];
PQC = Flatten[{Map[Subscript[Peq, #] &, Complement[GenIndex,
{SV}]],
Map[Subscript[Peq, #] &, LoadIndex],
Map[Subscript[Qeq, #] &, LoadIndex]}];
eq = (PQC - PQS0) /. Flatten[{SVA, KV}]; error = 1;
s1 = FindRoot[eq == 0, MapThread[{#1, #2} &, {var, varcal0}]];
PLC = Sum[
Subscript[Peq, i] /. Flatten[{KV, s1, SVA}], {i, 1,
Length[BUSdat1]}];
PGC = Map[Subscript[Peq, #] &, GenIndex] /. Flatten[{KV, s1, SVA}];
QGC = Map[Subscript[Qeq, #] &, GenIndex] /. Flatten[{KV, s1, SVA}];
PDC = Map[Chop[#] &,
Map[Subscript[Peq, #] &, LoadIndex] /. Flatten[{KV, s1, SVA}]];
QDC = Map[Chop[#] &,
Map[Subscript[Qeq, #] &, LoadIndex] /. Flatten[{KV, s1, SVA}]];];
In[18]:= apl = 1;
PG = apl {0.716, 1.63, 0.85};
PD = apl {0, 1.25, 0.9, 0, 1, 0};
QD = apl {0, 0.5, 0.3, 0, 0.35, 0};
QPratio = {0.5/1.25, 0.3/0.9, 1/0.35};
BUSdat = {{1, 2, PG[[1]], 0, 1.04, 0}, {2, 2, PG[[2]], 0, 1.025, 0},
{3, 2,
PG[[3]], 0, 1.025, 0}, {4, 1, -PD[[1]], -QD[[1]], 1, -0.1}, {5,
1, -PD[[2]], -QD[[2]], 1, 0}, {6, 1, -PD[[3]], -QD[[3]], 1, 0},
{7,
1, -PD[[4]], -QD[[4]], 1, 0}, {8, 1, -PD[[5]], -QD[[5]], 1, 0},
{9,
1, -PD[[6]], -QD[[6]], 1, 0}};
PowerFlow[BUSdat, linedata, 1]
In[25]:= s1
Out[25]= {Subscript[\[Theta], 2] -> 0.161967, Subscript[\[Theta], 3] -
> 0.0814153,
Subscript[\[Theta], 4] -> -0.0386902, Subscript[\[Theta], 5] ->
-0.0696178,
Subscript[\[Theta], 6] -> -0.0643572, Subscript[\[Theta], 7] ->
0.064921,
Subscript[\[Theta], 8] -> 0.0126979, Subscript[\[Theta], 9] ->
0.0343257,
Subscript[V, 4] -> 1.02579, Subscript[V, 5] -> 0.995631,
Subscript[V, 6] -> 1.01265, Subscript[V, 7] -> 1.02577,
Subscript[V, 8] -> 1.01588, Subscript[V, 9] -> 1.03235}
Machine Data
In[26]:= ws = 2 Pi 60;
MH = 2 {23.64, 6.4, 3.01}/ws;
Dp = {0, 0, 0};
Xd = {0.146, 0.8958, 1.3125};
Xdp = {0.0608, 0.1198, 0.1813};
Xq = {0.0969, 0.8645, 1.2578};
Xqp = {0.0969, 0.1969, 0.25};
Td0p = {8.96, 6, 5.89};
Tq0p = {0.31, 0.535, 0.6};
Rs = {0, 0, 0};
Exciter Data
In[36]:= KA = 2 {10, 10, 10};
TA = 2 {0.1, 0.1, 0.1};
KE = {1, 1, 1};
TE = {0.314, 0.314, 0.314};
KF = {0.063, 0.063, 0.063};
TF = {0.35, 0.35, 0.35};
SEfd = Map[0.0039 Exp[1.555 #] &,
Table[Subscript[Efdv, k], {k, 1, Length[GenIndex]}]];
Dynamic Equation on ith machine
In[43]:= ME[i_] := Module[{l, k}, k = Position[GenIndex, i][[1, 1]];
l = Chop[{( Subscript[\[Omega], i] - ws),
1/MH[[k]] (Subscript[Tmv, i] - Subscript[Eqpv, i]
Subscript[Iqv, i] -
Subscript[Edpv, i] Subscript[Idv,
i] - (Xqp[[k]] - Xdp[[k]]) Subscript[Idv, i] Subscript[Iqv,
i] -
Dp[[k]] ( Subscript[\[Omega], i] - ws)), ( -Subscript[Eqpv,
i] - (Xd[[k]] - Xdp[[k]]) Subscript[Idv, i] +
Subscript[Efdv, i])/
Td0p[[k]], ( -Subscript[Edpv, i] + (Xq[[k]] - Xqp[[k]])
Subscript[Iqv,
i])/Tq0p[[
k]], (-(KE[[k]] + SEfd[[k]]) Subscript[Efdv, i] +
Subscript[VRv, i])/
TE[[k]], 1/
TA[[k]] (-Subscript[VRv, i] +
KA[[k]] (Subscript[Vrefv, i] - Subscript[V, i] +
Subscript[Rfv, i] -
KF[[k]]/TF[[k]] Subscript[Efdv, i])), (-Subscript[Rfv,
i] + (KF[[k]]/TF[[k]]) Subscript[Efdv, i])/TF[[k]]},
10^-5];
Return[l];];
In[44]:= Arg1[x_] := If[x == 0, Return[0], Arg[x]];
Stator Equation on ith machine
In[45]:= SE[i_] := Module[{k, l}, k = Position[GenIndex, i][[1, 1]];
l = { Subscript[Edpv, i] -
Subscript[V, i] Sin[Subscript[\[Delta], i] - Subscript[\[Theta],
i]] -
Rs[[k]] Subscript[Idv, i] + Xqp[[k]] Subscript[Iqv, i],
Subscript[Eqpv, i] -
Subscript[V, i] Cos[Subscript[\[Delta], i] - Subscript[\[Theta],
i]] -
Rs[[k]] Subscript[Iqv, i] - Xdp[[k]] Subscript[Idv, i]};
Return[l]];
Network Equation on ith machine
In[46]:= NME[i_] :=
Module[{l},
l = {Subscript[Idv, i] Subscript[V, i]
Sin[Subscript[\[Delta], i] - Subscript[\[Theta], i]] +
Subscript[Iqv, i] Subscript[V, i]
Cos[Subscript[\[Delta], i] - Subscript[\[Theta], i]] - \!\(
\*UnderoverscriptBox[\(\[Sum]\), \(j = 1\), \(BUSN\)]
\*SubscriptBox[\(V\), \(i\)]\
\*SubscriptBox[\(V\), \(j\)]\ Abs[FY[\([\)\(i, j\)\(]\)]]\ Cos[
\*SubscriptBox[\(\[Theta]\), \(i\)] -
\*SubscriptBox[\(\[Theta]\), \(j\)] - Arg1[FY[\([\)\(i, j\)\(]\)]]]
\),
Subscript[Idv, i] Subscript[V, i]
Cos[Subscript[\[Delta], i] - Subscript[\[Theta], i]] -
Subscript[Iqv, i] Subscript[V, i]
Sin[Subscript[\[Delta], i] - Subscript[\[Theta], i]] - \!\(
\*UnderoverscriptBox[\(\[Sum]\), \(j = 1\), \(BUSN\)]
\*SubscriptBox[\(V\), \(i\)]\
\*SubscriptBox[\(V\), \(j\)]\ Abs[FY[\([\)\(i, j\)\(]\)]]\ Sin[
\*SubscriptBox[\(\[Theta]\), \(i\)] -
\*SubscriptBox[\(\[Theta]\), \(j\)] - Arg1[FY[\([\)\(i, j\)\(]\)]]]
\)};
Return[l];]
Network Equation on ith Load Bus
In[47]:= NLE[i_, BUSdat_] := Module[{l}, l = {BUSdat[[i, 3]] - \!\(
\*UnderoverscriptBox[\(\[Sum]\), \(k = 1\), \(BUSN\)]
\*SubscriptBox[\(V\), \(i\)]\
\*SubscriptBox[\(V\), \(k\)]\ Abs[FY[\([\)\(i, k\)\(]\)]]\ Cos[
\*SubscriptBox[\(\[Theta]\), \(i\)] -
\*SubscriptBox[\(\[Theta]\), \(k\)] - Arg1[FY[\([\)\(i, k\)\(]\)]]]
\),
BUSdat[[i, 4]] - \!\(
\*UnderoverscriptBox[\(\[Sum]\), \(k = 1\), \(BUSN\)]
\*SubscriptBox[\(V\), \(i\)]\
\*SubscriptBox[\(V\), \(k\)]\ Abs[FY[\([\)\(i, k\)\(]\)]]\ Sin[
\*SubscriptBox[\(\[Theta]\), \(i\)] -
\*SubscriptBox[\(\[Theta]\), \(k\)] - Arg1[FY[\([\)\(i, k\)\(]\)]]]
\)};
Return[l]];
Function calculating the initial values
In[48]:= TH = Flatten[{KV, SVA, s1}];
\[Omega]c =
Thread[Map[Subscript[\[Omega], #] &, GenIndex] ->
Table[ws, {i, 1, Length[GenIndex]}]];
IG = Map[((PGC[[Position[GenIndex, #][[1, 1]]]] -
I QGC[[Position[GenIndex, #][[1, 1]]]])/
Conjugate[Subscript[V, #] Exp[I Subscript[\[Theta], #]] /. TH])
&,
GenIndex];
Dlt = Map[
Arg[((Subscript[V, #]
Exp[I Subscript[\[Theta], #]] ) + (Rs[[
Position[GenIndex, #][[1, 1]]]] +
I Xq[[Position[GenIndex, #][[1, 1]]]])
IG[[Position[GenIndex, #][[1, 1]]]]) /. TH] &, GenIndex];
Dlt1 = 180/Pi Map[
Arg[((Subscript[V, #]
Exp[I Subscript[\[Theta], #]] ) + (Rs[[
Position[GenIndex, #][[1, 1]]]] +
I Xq[[Position[GenIndex, #][[1, 1]]]])
IG[[Position[GenIndex, #][[1, 1]]]]) /. TH] &, GenIndex];
Id = Map[Re[
IG[[Position[GenIndex, #][[1,
1]]]] Exp[-I (Dlt[[Position[GenIndex, #][[1, 1]]]] - Pi/2)]]
&,
GenIndex];
Iq = Map[Im[
IG[[Position[GenIndex, #][[1,
1]]]] Exp[-I (Dlt[[Position[GenIndex, #][[1, 1]]]] - Pi/2)]]
&,
GenIndex];
Vd = Map[Re[(Subscript[V, #]
Exp[I Subscript[\[Theta], #]] Exp[-I (Dlt[[
Position[GenIndex, #][[1, 1]]]] - Pi/2)]) /. TH] &,
GenIndex];
Vq = Map[Im[(Subscript[V, #]
Exp[I Subscript[\[Theta], #]] Exp[-I (Dlt[[
Position[GenIndex, #][[1, 1]]]] - Pi/2)]) /. TH] &,
GenIndex];
Edp = Table[(Xq[[k]] - Xqp[[k]]) Iq[[k]], {k, 1, Length[GenIndex]}];
Eqp = Map[Subscript[V, #]
Cos[Dlt[[Position[GenIndex, #][[1, 1]]]] - Subscript[\[Theta],
#]] +
Rs[[Position[GenIndex, #][[1, 1]]]] Iq[[
Position[GenIndex, #][[1, 1]]]] +
Xdp[[Position[GenIndex, #][[1, 1]]]] Id[[
Position[GenIndex, #][[1, 1]]]] &, GenIndex] /. TH;
Efd = Table[Eqp[[k]] + (Xd[[k]] - Xdp[[k]]) Id[[k]], {k, 1,
Length[Xdp]}];
SEfd2 = Map[0.0039 Exp[1.555 #] &, Efd];
VR = Table[(KE[[k]] + SEfd2[[k]]) Efd[[k]], {k, 1, Length[GenIndex]}];
Rf = Table[KF[[k]]/TF[[k]] Efd[[k]], {k, 1, Length[GenIndex]}];
Vref = Map[((Subscript[V, #] +
VR[[Position[GenIndex, #][[1, 1]]]]/
KA[[Position[GenIndex, #][[1, 1]]]]) /. TH) &, GenIndex];
Tms = Solve[
Map[(Subscript[Tm, Position[GenIndex, #][[1, 1]]] -
Eqp[[Position[GenIndex, #][[1, 1]]]] Iq[[
Position[GenIndex, #][[1,
1]]]] - (Xq[[Position[GenIndex, #][[1, 1]]]] -
Xdp[[Position[GenIndex, #][[1, 1]]]]) Id[[
Position[GenIndex, #][[1, 1]]]]
Iq[[Position[GenIndex, #][[1, 1]]]]) &, GenIndex] == 0,
Map[Subscript[Tm, Position[GenIndex, #][[1, 1]]] &, GenIndex]];
TM = Table[Tms[[1, i, 2]], {i, 1, Length[Tms[[1]]]}];
CV = Table[Subscript[V, i], {i, 1, BUSN}] /. TH;
CA = Table[Subscript[\[Theta], i], {i, 1, BUSN}] /. TH;
dltc = Thread[ Map[Subscript[\[Delta], #] &, GenIndex] -> Dlt];
dltc1 = Thread[ Map[Subscript[\[Delta], #] &, GenIndex] -> Dlt1];
Idc = Thread[Map[Subscript[Idv, #] &, GenIndex] -> Id];
Iqc = Thread[Map[Subscript[Iqv, #] &, GenIndex] -> Iq];
Vdc = Thread[Map[Subscript[Vdv, #] &, GenIndex] -> Vd];
Vqc = Thread[Map[Subscript[Vqv, #] &, GenIndex] -> Vq];
Eqpc = Thread[Map[Subscript[Eqpv, #] &, GenIndex] -> Eqp];
Edpc = Thread[Map[Subscript[Edpv, #] &, GenIndex] -> Edp];
Rfc = Thread[Map[Subscript[Rfv, #] &, GenIndex] -> Rf];
VRc = Thread[Map[Subscript[VRv, #] &, GenIndex] -> VR];
Efdc = Thread[Map[Subscript[Efdv, #] &, GenIndex] -> Efd];
Vrefc = Thread[Map[Subscript[Vrefv, #] &, GenIndex] -> Vref];
TMc = Thread[Map[Subscript[Tmv, #] &, GenIndex] -> TM];
dyc = Flatten[{dltc, Idc, Iqc, Edpc, Eqpc, Rfc, VRc, Efdc, Vrefc,
TMc, \[Omega]c}];
dyc1 = Flatten[{dltc1, Idc, Iqc, Vdc, Vqc, Edpc, Eqpc, Rfc, VRc, Efdc,
Vrefc,
TMc, \[Omega]c}];
Sort[dyc1, #1[[1, 2]] < #2[[1, 2]] &] // TableForm
Out[83]//TableForm= \!\(\*
TagBox[
TagBox[GridBox[{
{
RowBox[{
SubscriptBox["\[Omega]", "1"], "->",
RowBox[{"120", " ", "\[Pi]"}]}]},
{
RowBox[{
SubscriptBox["Tmv", "1"], "->", "0.7164102147448297`"}]},
{
RowBox[{
SubscriptBox["Vrefv", "1"], "->", "1.0952427423537974`"}]},
{
RowBox[{
SubscriptBox["Efdv", "1"], "->", "1.0821480400777919`"}]},
{
RowBox[{
SubscriptBox["VRv", "1"], "->", "1.1048548470759447`"}]},
{
RowBox[{
SubscriptBox["Rfv", "1"], "->", "0.19478664721400252`"}]},
{
RowBox[{
SubscriptBox["Eqpv", "1"], "->", "1.0563639528040119`"}]},
{
RowBox[{
SubscriptBox["Edpv", "1"], "->", "0.`"}]},
{
RowBox[{
SubscriptBox["Vqv", "1"], "->", "1.037964040758873`"}]},
{
RowBox[{
SubscriptBox["Vdv", "1"], "->", "0.06504344772160198`"}]},
{
RowBox[{
SubscriptBox["Iqv", "1"], "->", "0.6712430105428464`"}]},
{
RowBox[{
SubscriptBox["Idv", "1"], "->", "0.3026301323213611`"}]},
{
RowBox[{
SubscriptBox["\[Delta]", "1"], "->", "3.585720016443848`"}]},
{
RowBox[{
SubscriptBox["\[Omega]", "2"], "->",
RowBox[{"120", " ", "\[Pi]"}]}]},
{
RowBox[{
SubscriptBox["Tmv", "2"], "->", "1.63`"}]},
{
RowBox[{
SubscriptBox["Vrefv", "2"], "->", "1.1201038860909696`"}]},
{
RowBox[{
SubscriptBox["Efdv", "2"], "->", "1.7893233370769823`"}]},
{
RowBox[{
SubscriptBox["VRv", "2"], "->", "1.9020777218193912`"}]},
{
RowBox[{
SubscriptBox["Rfv", "2"], "->", "0.3220782006738568`"}]},
{
RowBox[{
SubscriptBox["Eqpv", "2"], "->", "0.7881690324629765`"}]},
{
RowBox[{
SubscriptBox["Edpv", "2"], "->", "0.6221979746111409`"}]},
{
RowBox[{
SubscriptBox["Vqv", "2"], "->", "0.6336093859516906`"}]},
{
RowBox[{
SubscriptBox["Vdv", "2"], "->", "0.8057072334501667`"}]},
{
RowBox[{
SubscriptBox["Iqv", "2"], "->", "0.9319921728746866`"}]},
{
RowBox[{
SubscriptBox["Idv", "2"], "->", "1.290147299760317`"}]},
{
RowBox[{
SubscriptBox["\[Delta]", "2"], "->", "61.098441521743766`"}]},
{
RowBox[{
SubscriptBox["\[Omega]", "3"], "->",
RowBox[{"120", " ", "\[Pi]"}]}]},
{
RowBox[{
SubscriptBox["Tmv", "3"], "->", "0.8500000000000001`"}]},
{
RowBox[{
SubscriptBox["Vrefv", "3"], "->", "1.097573933543132`"}]},
{
RowBox[{
SubscriptBox["Efdv", "3"], "->", "1.4029943029585485`"}]},
{
RowBox[{
SubscriptBox["VRv", "3"], "->", "1.451478670862643`"}]},
{
RowBox[{
SubscriptBox["Rfv", "3"], "->", "0.2525389745325387`"}]},
{
RowBox[{
SubscriptBox["Eqpv", "3"], "->", "0.7678611247614691`"}]},
{
RowBox[{
SubscriptBox["Edpv", "3"], "->", "0.6242375949756084`"}]},
{
RowBox[{
SubscriptBox["Vqv", "3"], "->", "0.666066883948942`"}]},
{
RowBox[{
SubscriptBox["Vdv", "3"], "->", "0.779089151578012`"}]},
{
RowBox[{
SubscriptBox["Iqv", "3"], "->", "0.6194062264096134`"}]},
{
RowBox[{
SubscriptBox["Idv", "3"], "->", "0.561468509721605`"}]},
{
RowBox[{
SubscriptBox["\[Delta]", "3"], "->", "54.136617534872826`"}]}
},
GridBoxAlignment->{
"Columns" -> {{Left}}, "ColumnsIndexed" -> {}, "Rows" ->
{{Baseline}},
"RowsIndexed" -> {}},
GridBoxSpacings->{"Columns" -> {
Offset[0.27999999999999997`], {
Offset[0.5599999999999999]},
Offset[0.27999999999999997`]}, "ColumnsIndexed" -> {}, "Rows" -> {
Offset[0.2], {
Offset[0.4]},
Offset[0.2]}, "RowsIndexed" -> {}}],
Column],
Function[BoxForm`e$,
TableForm[BoxForm`e$]]]\)
Variable Definition
In[84]:= dyvar = Flatten[
Table[{Subscript[\[Delta], i], Subscript[\[Omega], i],
Subscript[Eqpv, i],
Subscript[Edpv, i], Subscript[Efdv, i], Subscript[VRv, i],
Subscript[Rfv,
i]}, {i, 1, Length[GenIndex]}]];
Idqvar = Flatten[Map[{Subscript[Idv, #], Subscript[Iqv, #]} &,
GenIndex]];
Vg = Flatten[Map[{Subscript[\[Theta], #], Subscript[V, #]} &,
GenIndex]];
Vl = Flatten[Map[{Subscript[\[Theta], #], Subscript[V, #]} &,
LoadIndex]];
Uc = Flatten[Map[{Subscript[Tmv, #], Subscript[Vrefv, #]} &,
GenIndex]];
Gvar = Flatten[Map[{Subscript[\[Theta], #]} &, Complement[GenIndex,
{SV}]]];
The Whole Equations
In[90]:= Off[General::"spell"]
WDE = Chop[Flatten[Map[ME[#] &, GenIndex]]];(*the whole dynamic
equations*)
WSE = Chop[Flatten[Map[SE[#] &, GenIndex]]];(*the whole stator
equations*);
WNLE = Chop[
Flatten[Map[NLE[#, BUSdat] &,
LoadIndex]]];(*the whole network load equations*)
WNME = Chop[
Flatten[Map[NME[#] &, GenIndex]]];(*the whole network machine
equations*)
ic = Flatten[{dyc, TH}];
te2 = Flatten[{WDE, WSE, WNLE, WNME}];
Chop[te2 /. ic]
Out[97]= {0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0,
0, 0, 0, 0, 0, \
0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0}
Eigenvalue sensitivity
Matrix manipulation 1
In[98]:= A15 = Chop[Outer[D, WDE, dyvar]];
B15 = Chop[Outer[D, WDE, Idqvar]];
B25 = Chop[Outer[D, WDE, Vg]];
C15 = Chop[Outer[D, WSE, dyvar]];
D15 = Chop[Outer[D, WSE, Idqvar]];
D25 = Chop[Outer[D, WSE, Vg]];
C25 = Chop[Outer[D, WNME, dyvar]];
D35 = Chop[Outer[D, WNME, Idqvar]];
D45 = Chop[Outer[D, WNME, Vg]];
D55 = Chop[Outer[D, WNME, Vl]];
D65 = Chop[Outer[D, WNLE, Vg]];
D75 = Chop[Outer[D, WNLE, Vl]];
Matrix manipulation 2
In[110]:= (*Model Reduction*)
Ap5 = (A15 - B15.Inverse[D15].C15);
Bp5 = (B25 - B15.Inverse[D15].D25);
K15 = D45 - D35.Inverse[D15].D25;
K25 = C25 - D35.Inverse[D15].C15;
Matrix manipulation 3
(*My own way*)
temp1 = K15 - D55.Inverse[D75].D65;
Timing[Asys51pp = Ap5 - Bp5.Inverse[temp1].K25];
Asyspssonly = Asys51pp;