From 68c835ca1e79337727c3104718b1e4285595f3ce Mon Sep 17 00:00:00 2001 From: Yang Liu Date: Tue, 25 Apr 2023 03:22:42 -0400 Subject: [PATCH] enable emt in runPowerSAS.m function --- EMT/calEMT.m | 6 +- EMT/plotEMTResults.m | 362 ++++++++++++++++---------------- example/ex_emt.m | 7 +- example/output/~default_log.txt | 6 + runPowerSAS.m | 2 + 5 files changed, 192 insertions(+), 191 deletions(-) create mode 100644 example/output/~default_log.txt diff --git a/EMT/calEMT.m b/EMT/calEMT.m index 74d45b9..270afb6 100644 --- a/EMT/calEMT.m +++ b/EMT/calEMT.m @@ -1,5 +1,5 @@ -function [t_final, sol_final] = calEMT(testSystem) +function [res] = calEMT(testSystem) %% [orderOfDT, timeStepOfDT] = setSASParameters(); @@ -49,6 +49,6 @@ x_ini_post = x_on(end,:)'; t_post = t_dt; x_post = x_dt'; -t_final = [t_pre; t_on(2:end); t_post(2:end)]; -sol_final = [x_pre; x_on(2:end,:); x_post(2:end,:)]; +res.t_final = [t_pre; t_on(2:end); t_post(2:end)]; +res.sol_final = [x_pre; x_on(2:end,:); x_post(2:end,:)]; return \ No newline at end of file diff --git a/EMT/plotEMTResults.m b/EMT/plotEMTResults.m index 13b77eb..aa29153 100644 --- a/EMT/plotEMTResults.m +++ b/EMT/plotEMTResults.m @@ -1,97 +1,93 @@ +function plotEMTResults(res) close all - -% t_final = [t_pre; t_on(2:end); t_post(2:end)]; -% sol_final = [x_pre; x_on(2:end,:); x_post(2:end,:)]; -% % t_final = [t_pre; t_on(2:end)]; -% % sol_final = [x_pre; x_on(2:end,:)]; -% % t_final = [t_pre]; -% % sol_final = [x_pre]; - -%% relative rotor angles -figure(1) -% EMT result -plot(t_final, (sol_final(:,1)-sol_final(:,1) )*180/pi,'r--','LineWidth',3) -hold on -plot(t_final, (sol_final(:,2)-sol_final(:,1) )*180/pi,'g--','LineWidth',3) -hold on -plot(t_final, (sol_final(:,3)-sol_final(:,1) )*180/pi,'b--','LineWidth',3) -hold on -plot(t_final, (sol_final(:,4)-sol_final(:,1) )*180/pi,'b--','LineWidth',3) -hold on -% ylim([-60,10]) -xlabel('Time (s)') -ylabel('Relative rotor angles (degree)') -set(gca, 'Fontname', 'Times New Roman', 'Fontsize', 12); -title('Relative rotor angles') -legend('EMT G1','EMT G2','EMT G3','EMT G4') - -%% frequency deviation -figure(2) -% EMT result -plot(t_final, sol_final(:,9),'r--','LineWidth',3) -hold on -plot(t_final, sol_final(:,10),'g--','LineWidth',3) -hold on -plot(t_final, sol_final(:,11),'b--','LineWidth',3) -hold on -plot(t_final, sol_final(:,12),'m--','LineWidth',3) -hold on -% ylim([-60,10]) -% xlim([0,0.25]) -xlabel('Time (s)') -ylabel('Frequency deviation (p.u.)') -set(gca, 'Fontname', 'Times New Roman', 'Fontsize', 12); -title('Frequency deviation') -legend('EMT G1','EMT G2','EMT G3','EMT G4') - -%% Exciter field voltage -figure(3) -% EMT result -plot(t_final, sol_final(:,33),'r--','LineWidth',3) -hold on -plot(t_final, sol_final(:,34),'g--','LineWidth',3) -hold on -plot(t_final, sol_final(:,35),'b--','LineWidth',3) -hold on -plot(t_final, sol_final(:,36),'m--','LineWidth',3) -xlabel('Time (s)') -xlim([0,0.4]) -ylabel('Exciter field voltage (p.u.)') -set(gca, 'Fontname', 'Times New Roman', 'Fontsize', 12); -title('Exciter field voltage') -legend('EMT G1','EMT G2','EMT G3','EMT G4') - -%% P1 -figure(4) -% EMT result -plot(t_final, sol_final(:,37),'r--','LineWidth',3) -hold on -plot(t_final, sol_final(:,38),'g--','LineWidth',3) -hold on -plot(t_final, sol_final(:,39),'b--','LineWidth',3) -hold on -plot(t_final, sol_final(:,40),'m--','LineWidth',3) - -xlabel('Time (s)') -ylabel('P1(p.u.)') -set(gca, 'Fontname', 'Times New Roman', 'Fontsize', 12); -title('P1') -legend('EMT G1','EMT G2','EMT G3','EMT G4') - -%% Voltage at bus 7 -figure(5) -% EMT result -plot(t_final, sol_final(:,51),'r','LineWidth',1) -hold on -plot(t_final, sol_final(:,52),'g','LineWidth',1) -hold on -plot(t_final, sol_final(:,53),'b','LineWidth',1) - -xlabel('Time (s)') -ylabel('Voltage at bus 7(p.u.)') -set(gca, 'Fontname', 'Times New Roman', 'Fontsize', 12); -title('Voltage at bus 7') -legend('EMT phase A','EMT phase B','EMT phase C') +t_final = res.t_final; +sol_final = res.sol_final; + +% %% relative rotor angles +% figure(1) +% % EMT result +% plot(t_final, (sol_final(:,1)-sol_final(:,1) )*180/pi,'r--','LineWidth',3) +% hold on +% plot(t_final, (sol_final(:,2)-sol_final(:,1) )*180/pi,'g--','LineWidth',3) +% hold on +% plot(t_final, (sol_final(:,3)-sol_final(:,1) )*180/pi,'b--','LineWidth',3) +% hold on +% plot(t_final, (sol_final(:,4)-sol_final(:,1) )*180/pi,'b--','LineWidth',3) +% hold on +% % ylim([-60,10]) +% xlabel('Time (s)') +% ylabel('Relative rotor angles (degree)') +% set(gca, 'Fontname', 'Times New Roman', 'Fontsize', 12); +% title('Relative rotor angles') +% legend('EMT G1','EMT G2','EMT G3','EMT G4') + +% %% frequency deviation +% figure(2) +% % EMT result +% plot(t_final, sol_final(:,9),'r--','LineWidth',3) +% hold on +% plot(t_final, sol_final(:,10),'g--','LineWidth',3) +% hold on +% plot(t_final, sol_final(:,11),'b--','LineWidth',3) +% hold on +% plot(t_final, sol_final(:,12),'m--','LineWidth',3) +% hold on +% % ylim([-60,10]) +% % xlim([0,0.25]) +% xlabel('Time (s)') +% ylabel('Frequency deviation (p.u.)') +% set(gca, 'Fontname', 'Times New Roman', 'Fontsize', 12); +% title('Frequency deviation') +% legend('EMT G1','EMT G2','EMT G3','EMT G4') + +% %% Exciter field voltage +% figure(3) +% % EMT result +% plot(t_final, sol_final(:,33),'r--','LineWidth',3) +% hold on +% plot(t_final, sol_final(:,34),'g--','LineWidth',3) +% hold on +% plot(t_final, sol_final(:,35),'b--','LineWidth',3) +% hold on +% plot(t_final, sol_final(:,36),'m--','LineWidth',3) +% xlabel('Time (s)') +% xlim([0,0.4]) +% ylabel('Exciter field voltage (p.u.)') +% set(gca, 'Fontname', 'Times New Roman', 'Fontsize', 12); +% title('Exciter field voltage') +% legend('EMT G1','EMT G2','EMT G3','EMT G4') + +% %% P1 +% figure(4) +% % EMT result +% plot(t_final, sol_final(:,37),'r--','LineWidth',3) +% hold on +% plot(t_final, sol_final(:,38),'g--','LineWidth',3) +% hold on +% plot(t_final, sol_final(:,39),'b--','LineWidth',3) +% hold on +% plot(t_final, sol_final(:,40),'m--','LineWidth',3) +% +% xlabel('Time (s)') +% ylabel('P1(p.u.)') +% set(gca, 'Fontname', 'Times New Roman', 'Fontsize', 12); +% title('P1') +% legend('EMT G1','EMT G2','EMT G3','EMT G4') + +% %% Voltage at bus 7 +% figure(5) +% % EMT result +% plot(t_final, sol_final(:,51),'r','LineWidth',1) +% hold on +% plot(t_final, sol_final(:,52),'g','LineWidth',1) +% hold on +% plot(t_final, sol_final(:,53),'b','LineWidth',1) +% +% xlabel('Time (s)') +% ylabel('Voltage at bus 7(p.u.)') +% set(gca, 'Fontname', 'Times New Roman', 'Fontsize', 12); +% title('Voltage at bus 7') +% legend('EMT phase A','EMT phase B','EMT phase C') %% Voltage at bus 8 figure(6) @@ -108,50 +104,50 @@ set(gca, 'Fontname', 'Times New Roman', 'Fontsize', 12); title('Voltage at bus 8') legend('EMT phase A','EMT phase B','EMT phase C') -%% Voltage at bus 10 -figure(7) -% EMT result -plot(t_final, sol_final(:,60),'r','LineWidth',1) -hold on -plot(t_final, sol_final(:,61),'g','LineWidth',1) -hold on -plot(t_final, sol_final(:,62),'b','LineWidth',1) - -xlabel('Time (s)') -ylabel('Voltage at bus 10(p.u.)') -set(gca, 'Fontname', 'Times New Roman', 'Fontsize', 12); -title('Voltage at bus 10') -legend('EMT phase A','EMT phase B','EMT phase C') - -%% Voltage at bus 6 -figure(8) -% EMT result -plot(t_final, sol_final(:,48),'r','LineWidth',1) -hold on -plot(t_final, sol_final(:,49),'g','LineWidth',1) -hold on -plot(t_final, sol_final(:,50),'b','LineWidth',1) - -xlabel('Time (s)') -ylabel('Voltage at bus 6(p.u.)') -set(gca, 'Fontname', 'Times New Roman', 'Fontsize', 12); -title('Voltage at bus 6') -legend('EMT phase A','EMT phase B','EMT phase C') - -%% current on branch 5-6 -figure(9) -% EMT result -plot(t_final, sol_final(:,78),'r','LineWidth',1) -hold on -plot(t_final, sol_final(:,79),'g','LineWidth',1) -hold on -plot(t_final, sol_final(:,80),'b','LineWidth',1) - -xlabel('Time (s)') -ylabel('Current on branch 5-6(p.u.)') -set(gca, 'Fontname', 'Times New Roman', 'Fontsize', 12); -title('Current on branch 5-6') -legend('EMT phase A','EMT phase B','EMT phase C') +% %% Voltage at bus 10 +% figure(7) +% % EMT result +% plot(t_final, sol_final(:,60),'r','LineWidth',1) +% hold on +% plot(t_final, sol_final(:,61),'g','LineWidth',1) +% hold on +% plot(t_final, sol_final(:,62),'b','LineWidth',1) +% +% xlabel('Time (s)') +% ylabel('Voltage at bus 10(p.u.)') +% set(gca, 'Fontname', 'Times New Roman', 'Fontsize', 12); +% title('Voltage at bus 10') +% legend('EMT phase A','EMT phase B','EMT phase C') + +% %% Voltage at bus 6 +% figure(8) +% % EMT result +% plot(t_final, sol_final(:,48),'r','LineWidth',1) +% hold on +% plot(t_final, sol_final(:,49),'g','LineWidth',1) +% hold on +% plot(t_final, sol_final(:,50),'b','LineWidth',1) +% +% xlabel('Time (s)') +% ylabel('Voltage at bus 6(p.u.)') +% set(gca, 'Fontname', 'Times New Roman', 'Fontsize', 12); +% title('Voltage at bus 6') +% legend('EMT phase A','EMT phase B','EMT phase C') + +% %% current on branch 5-6 +% figure(9) +% % EMT result +% plot(t_final, sol_final(:,78),'r','LineWidth',1) +% hold on +% plot(t_final, sol_final(:,79),'g','LineWidth',1) +% hold on +% plot(t_final, sol_final(:,80),'b','LineWidth',1) +% +% xlabel('Time (s)') +% ylabel('Current on branch 5-6(p.u.)') +% set(gca, 'Fontname', 'Times New Roman', 'Fontsize', 12); +% title('Current on branch 5-6') +% legend('EMT phase A','EMT phase B','EMT phase C') %% current on branch 6-7 figure(10) @@ -168,49 +164,49 @@ set(gca, 'Fontname', 'Times New Roman', 'Fontsize', 12); title('Current on branch 6-7') legend('EMT phase A','EMT phase B','EMT phase C') -%% current on branch 7-8 -figure(11) -% EMT result -plot(t_final, sol_final(:,84),'r','LineWidth',1) -hold on -plot(t_final, sol_final(:,85),'g','LineWidth',1) -hold on -plot(t_final, sol_final(:,86),'b','LineWidth',1) - -xlabel('Time (s)') -ylabel('Current on branch 7-8(p.u.)') -set(gca, 'Fontname', 'Times New Roman', 'Fontsize', 12); -title('Current on branch 7-8') -legend('EMT phase A','EMT phase B','EMT phase C') - -%% current on branch 9-10 -figure(12) -% EMT result -plot(t_final, sol_final(:,96),'r','LineWidth',1) -hold on -plot(t_final, sol_final(:,97),'g','LineWidth',1) -hold on -plot(t_final, sol_final(:,98),'b','LineWidth',1) - -xlabel('Time (s)') -ylabel('Current on branch 9-10(p.u.)') -set(gca, 'Fontname', 'Times New Roman', 'Fontsize', 12); -title('Current on branch 9-10') -legend('EMT phase A','EMT phase B','EMT phase C') - -%% Voltage at bus 6 -figure(13) -% EMT result -plot(t_final, sol_final(:,108),'r','LineWidth',1) -hold on -plot(t_final, sol_final(:,109),'g','LineWidth',1) -hold on -plot(t_final, sol_final(:,110),'b','LineWidth',1) -hold on -plot(t_final, sol_final(:,111),'b','LineWidth',1) - -xlabel('Time (s)') -ylabel('Voltage at bus 6(p.u.)') -set(gca, 'Fontname', 'Times New Roman', 'Fontsize', 12); -title('Voltage at bus 6') -legend('Gen1','Gen2','Gen3','Gen4') +% %% current on branch 7-8 +% figure(11) +% % EMT result +% plot(t_final, sol_final(:,84),'r','LineWidth',1) +% hold on +% plot(t_final, sol_final(:,85),'g','LineWidth',1) +% hold on +% plot(t_final, sol_final(:,86),'b','LineWidth',1) +% +% xlabel('Time (s)') +% ylabel('Current on branch 7-8(p.u.)') +% set(gca, 'Fontname', 'Times New Roman', 'Fontsize', 12); +% title('Current on branch 7-8') +% legend('EMT phase A','EMT phase B','EMT phase C') + +% %% current on branch 9-10 +% figure(12) +% % EMT result +% plot(t_final, sol_final(:,96),'r','LineWidth',1) +% hold on +% plot(t_final, sol_final(:,97),'g','LineWidth',1) +% hold on +% plot(t_final, sol_final(:,98),'b','LineWidth',1) +% +% xlabel('Time (s)') +% ylabel('Current on branch 9-10(p.u.)') +% set(gca, 'Fontname', 'Times New Roman', 'Fontsize', 12); +% title('Current on branch 9-10') +% legend('EMT phase A','EMT phase B','EMT phase C') + +% %% Voltage at bus 6 +% figure(13) +% % EMT result +% plot(t_final, sol_final(:,108),'r','LineWidth',1) +% hold on +% plot(t_final, sol_final(:,109),'g','LineWidth',1) +% hold on +% plot(t_final, sol_final(:,110),'b','LineWidth',1) +% hold on +% plot(t_final, sol_final(:,111),'b','LineWidth',1) +% +% xlabel('Time (s)') +% ylabel('Voltage at bus 6(p.u.)') +% set(gca, 'Fontname', 'Times New Roman', 'Fontsize', 12); +% title('Voltage at bus 6') +% legend('Gen1','Gen2','Gen3','Gen4') diff --git a/example/ex_emt.m b/example/ex_emt.m index d64f6d9..1b50696 100644 --- a/example/ex_emt.m +++ b/example/ex_emt.m @@ -1,10 +1,7 @@ - clc; clear; close; addpath('EMT') -testSystem = 'twoarea'; +res = runPowerSAS('emt','twoarea'); -[t_final, sol_final] = calEMT(testSystem); - -plotEMTResults(); +plotEMTResults(res); diff --git a/example/output/~default_log.txt b/example/output/~default_log.txt new file mode 100644 index 0000000..2da6ae8 --- /dev/null +++ b/example/output/~default_log.txt @@ -0,0 +1,6 @@ + +03:20:33[INFO]The session ID is 20230425T032033017. +03:20:33[INFO]CASE twoarea. Timestamp=20230425T032033017. +03:20:33[INFO]Finished loading data file. +03:20:41[INFO] + diff --git a/runPowerSAS.m b/runPowerSAS.m index 80a88ef..a39652d 100644 --- a/runPowerSAS.m +++ b/runPowerSAS.m @@ -99,6 +99,8 @@ end %% main body if strcmp(simType,'pf') res=calPF(SysData,options,caseName); +elseif strcmp(simType,'emt') + res=calEMT(data); elseif strcmp(simType,'cpf') res=calCPF(SysData,options,caseName,varargin{:}); elseif strcmp(simType,'tsa')