]> 2012-02-07T15:02:48+00:00 NSW Electricity Network '330 kV' 1 '1.30 AUD 2009/watt' 0.945 0.96 '0.906t/MWh' 1 '100 MVA' '9999 MW' 1.1 '9999 MW' 0.9 %Graphical Controls %P_N = 4; %Number of plots (ONE OR FOUR) %SCENARIO_NUMBERS = [ 1, 3, 4, 5]; %DEFINE SCENARIOS TO COMPARE(between 0 and 5) IN MULTIPLOT. For single plot, first scenario is %V_P = [-33.75,150.5,8]; %Define map viewpoint (lat, long, zoom number). Zoom number MUST be between 6 & displayed. L_P = [25,100]; % Legend Position ([0,0] = TOP LEFT) for PLOTS = 1:P_N; mpc.version = '2'; %Get the NSW network defaults baseMVA %MachineBaseMva ?MaximumRealPowerOutput ?MinimumRealPowerOutput ?MaximumVoltageMagnitude ?MinimumVoltageMagnitude ?BusArea ?VoltageAngle ?BaseVoltage ?ShuntConductance nswDefaults = sparql("SELECT ?MachineBaseMva ?MaximumRealPowerOutput ?MinimumRealPowerOutput ?MaximumVoltageMagnitude ?MinimumVoltageMagnitude ?BusArea ?VoltageAngle ?BaseVoltage ?ShuntConductance WHERE{ a:NSW_Electricity_Network prop:MachineBaseMva ?MachineBaseMva. a:NSW_Electricity_Network prop:MaximumRealPowerOutput ?MaximumRealPowerOutput. a:NSW_Electricity_Network prop:MinimumRealPowerOutput ?MinimumRealPowerOutput. a:NSW_Electricity_Network prop:MaximumVoltageMagnitude ?MaximumVoltageMagnitude. a:NSW_Electricity_Network prop:MinimumVoltageMagnitude ?MinimumVoltageMagnitude. a:NSW_Electricity_Network prop:BusArea ?BusArea. a:NSW_Electricity_Network prop:VoltageAngle ?VoltageAngle. a:NSW_Electricity_Network prop:BaseVoltage ?BaseVoltage. a:NSW_Electricity_Network prop:ShuntConductance ?ShuntConductance. }"); %Parse the baseMVA mpc.baseMVA = sscanf (char(nswDefaults(1,1)), "%f MVA"); %?bus ?id ?BusType ?realPowerDemand ?reactivPowerDemand ?shuntSusceptance ?VoltageMagnitude ?RealPowerOutput ?MaximumReactivePowerOutput ?MinimumReactivePowerOutput ?VoltageMagnitudeSetpoint ?GeoLocation ?ShortName bus = sparql("select ?bus ?id ?BusType ?realPowerDemand ?reactivPowerDemand ?shuntSusceptance ?VoltageMagnitude ?RealPowerOutput ?MaximumReactivePowerOutput ?MinimumReactivePowerOutput ?VoltageMagnitudeSetpoint ?GeoLocation ?ShortName WHERE{ ?bus rdf:type cat:Electrical_Bus. OPTIONAL{?bus prop:BusIdentifier ?id.} OPTIONAL{?bus prop:BusType ?BusType.} OPTIONAL{?bus prop:RealPowerDemand ?realPowerDemand.} OPTIONAL{?bus prop:ReactivPowerDemand ?reactivPowerDemand.} OPTIONAL{?bus prop:ShuntSusceptance ?shuntSusceptance.} OPTIONAL{?bus prop:VoltageMagnitude ?VoltageMagnitude.} OPTIONAL{?bus prop:RealPowerOutput ?RealPowerOutput.} OPTIONAL{?bus prop:MaximumReactivePowerOutput ?MaximumReactivePowerOutput.} OPTIONAL{?bus prop:MinimumReactivePowerOutput ?MinimumReactivePowerOutput.} OPTIONAL{?bus prop:VoltageMagnitudeSetpoint ?VoltageMagnitudeSetpoint.} OPTIONAL{?bus prop:GeoLocation ?GeoLocation.} OPTIONAL{?bus prop:ShortName ?ShortName.} } ORDER BY ?bus"); [m,n] =size(bus); B_n = n; busCtr = 0; genCtr = 0; %----------------------------------------------------------------------------------------------------------------% %Build up the mpc.bus tablefrom bus. Take defaults from the nswDefaults variable [Bus Data (mpc.bus) matrix definition ] for i = 1:n BUS_I(i) = cell2mat(bus(2,i)); %BUS_I (column 1) Bus number BUS_TYPE(i) = cell2mat(bus(3,i)); %BUS_TYPE (colomn 2) Bus type (1 = PQ, 2 = PV, 3 = ref, 4 = isolated) PD(i) =sscanf (char(bus(4,i)), "%d MW"); %PD (column 3) Real Power Demand QD(i) =sscanf (char(bus(5,i)), "%d MVA"); %QD (column 4) Reactive Power Demand GS(i) = sscanf (char(nswDefaults(9,1)), "%f MW"); %GS (Column 5) Shunt Conductance if (length (char(bus(6,i)))<4), BS(i) = 0; %BS (Column 6) Shunt Susceptance else BS(i) = sscanf (char(bus(6,i)), "%f MVA"); endif; BUS_AREA(i) = cell2mat(nswDefaults(6,1)); %BUS_AREA (Column 7) Area Number if (length (char(bus(7,i)))<1), VM(i) = 1; else VM(i) = cell2mat(bus(7,i)); %VM (Column 8) Volatge Magnitude endif; VA(i) = cell2mat(nswDefaults(7,1)); %VA (Column 9) Volatge Angle BASE_KV(i) = sscanf (char(nswDefaults(8,1)), "%f kV"); %BASE_KV (Column 10) Base Volatge (kV) ZONE(i) = 1; %ZONE (Column 11) Loss Zone VMAX(i) = cell2mat(nswDefaults(4,1)); %VMAX (Coulmn 12) Maximum voltage magnitude VMIN(i) = cell2mat(nswDefaults(5,1)); %VMIN (Coulmn 13) Minimum voltage magnitude endfor mpc.bus = [BUS_I', BUS_TYPE', PD',QD',GS',BS',BUS_AREA',VM',VA', BASE_KV',ZONE',VMAX',VMIN']; %----------------------------------------------------------------------------------------------------------------% %Build up mpc.gen table from bus. Take defaults from the nswDefaults variable for j = 1:n if cell2mat(bus(3,j)) ==1, GEN_BUS(j)=0; %GEN_BUS (Column 1) Bus number else GEN_BUS(j)=cell2mat(bus(2,j)); endif if (length (char(bus(8,j)))<4), PG(j) = 0; %PG (Column 2) Real Power Output else PG(j) = sscanf (char(bus(8,j)), "%f MW"); endif; QG(j) = 0; %QG (Column 3) Reactive Power Output !HARD CODED if (length (char(bus(9,j)))<4), QMAX(j) = 0; %QMAX (Column 4) Maximum Reactive Power Output else QMAX(j) = sscanf (char(bus(9,j)), "%f MVA"); endif; if (length (char(bus(10,j)))<4), QMIN(j) = 0; %QMIN (Column 5) Minimum Reactive Power Output else QMIN(j) = sscanf (char(bus(10,j)), "%f MVA"); endif; if (length (char(bus(11,j)))<1), VG(j) = 0; %VG (Column 6) Voltage Magnitude else VG(j) = cell2mat(bus(11,j)); endif; MBASE(j) = sscanf (char(nswDefaults(1,1)), "%f MVA"); %MBASE (Column 7) Total MVA base of machine !HARDCODED GEN_STATUS(j) = 1; %GEN_STATUS (Column 8) Machine Status !HARDCODED PMAX(j) = 9999; %PMAX (Column 9) Maximum reaal power output !HARDCODED PMIN(j) = 9999; %PMIN (Column 10) Minumum real power output !HARDCODED endfor mpc.gen = [GEN_BUS', PG',QG', QMAX', QMIN', VG',MBASE',GEN_STATUS',PMAX',PMIN']; a = mpc.gen; [s, i] = sort (a (:, 1)); mpc.gen = a (i, :); len = length (mpc.gen)-length (mpc.gen(:,1)(mpc.gen(:,1)~=0)); %Delete zero leading rows mpc.gen(1:len,:)=[]; %----------------------------------------------------------------------------------------------------------------% %Get the line data %?line ?ShortName ?FromBus ?FromId ?FromGeoLocation ?ToBus ?ToId ?ToGeoLocation ?Resistance ?Reactance ?TotalLineChargingSusceptance branch = sparql("select ?line ?ShortName ?FromBus ?FromId ?FromGeoLocation ?ToBus ?ToId ?ToGeoLocation ?Resistance ?Reactance ?TotalLineChargingSusceptance WHERE{ ?line rdf:type cat:Electrical_Transmission. ?line prop:ShortName ?ShortName. ?line prop:FromBus ?FromBus. ?FromBus prop:BusIdentifier ?FromId. ?FromBus prop:GeoLocation ?FromGeoLocation. ?line prop:ToBus ?ToBus. ?ToBus prop:BusIdentifier ?ToId. ?ToBus prop:GeoLocation ?ToGeoLocation. ?line prop:Resistance ?Resistance. ?line prop:Reactance ?Reactance. ?line prop:TotalLineChargingSusceptance ?TotalLineChargingSusceptance . } ORDER BY ?line"); %build up the mpc.branch matrix [m,n] =size(branch); for i = 1:n FBUS(i) = cell2mat(branch(4,i)); %F_Bus From Bus (column 1) TBUS(i) = cell2mat(branch(7,i)); %T_BUS To Bus (column 2) BR_R(i) = cell2mat(branch(9,i)); %BR_R Resistance (column 3) BR_X(i) = cell2mat(branch(10,i)); %BR_X Reactance (column 4) BR_B(i) = cell2mat(branch(11,i)); %BR_B Total line charging suseptance (column 5) RATE_A(i) = 1000; %RATE_A MVA rating A (column 6) !HARDCODED RATE_B(i) = 1000; %RATE_B MVA rating B (column 7) !HARDCODED RATE_C(i) = 1000; %RATE_C MVA rating C (column 8) !HARDCODED TAP(i) = 0; %TAP Transformer off Nominal turns ratio (column 9) !HARDCODED SHIFT(i) = 0; %SHIFT Transformer phase shift angles (column 10) !HARDCODED BR_STATUS(i) = 1; %BR_STATUS Initial branch status (column 11) !HARDCODED ANGMIN(i) = -360; %ANGMIN Minimum angle difference (column 12) !HARDCODED ANGMAX(i) = 360; %ANDMAX Maximum angle difference (column 13) !HARDCODED endfor mpc.branch = [FBUS',TBUS',BR_R',BR_X',BR_B', RATE_A', RATE_B', RATE_C', TAP',SHIFT',BR_STATUS',ANGMIN',ANGMAX']; mpc; %----------------------------------------------------------------------------------------------------------------% %Get the scenario data %----------------------------------------------------------------------------------------------------------------% %Evaluate the relevant case SCENARIO = [ 1, 0.6, 1.05, 1.05, 1.05, 1.05; %SCENARIO MATRIX 1, 0.6, 1.05, 1.05, 1.05, 1.05; 1, 1, 1, 0, 1, 1; 1, 1, 1, 1, 0, 1; 1, 1, 1, 1, 1, 0 ]; SN = SCENARIO_NUMBERS(1,PLOTS)+1; % Scenario 0 - No modification % Scenario 1 - Base case (light load) % Scenario 2 - Heavy load % Scenario 3 - Heavy loading plus outage of the lines 12-13 and 12-16 (Vales - Sydney North Vales - Munmorah) % Scenario 4 - Heavy loading plus outage of generator #6 (Liddell) % Scenario 5 - Heavy loading plus outage of the lines 9-19 and 9-20 (Bayswater - Regentville and Bayswater - Mt Piper) mpc.bus(:,[3:4]) = SCENARIO(1,SN)*mpc.bus(:,[3:4]); mpc.gen(:,[2]) = SCENARIO(2,SN)*mpc.gen(:,[2]); mpc.branch([18,19],[11]) = SCENARIO(3,SN)*mpc.branch([18,19],[11]); mpc.gen([2],[8]) = SCENARIO(4,SN)*mpc.gen([2],[8]); mpc.branch([27,48],[11]) = SCENARIO(5,SN)*mpc.branch([27,48],[11]); diary off; results = runpf(mpc, mpoption( "OUT_ALL",0)); diary on; %----------------------------------------------------------------------------------------------------------------% %Apparent Power Calculation PP_IN = diag(results.branch(:,14)*(results.branch(:,14))'); %Test Apparent power in is under limit QQ_IN = diag(results.branch(:,15)*(results.branch(:,15))'); APPARENT_POWER_IN = sqrt (PP_IN+QQ_IN); PP_OUT = diag(results.branch(:,16)*(results.branch(:,16))'); %Test Apparent power out is under limit QQ_OUT = diag(results.branch(:,17)*(results.branch(:,17))'); APPARENT_POWER_OUT = sqrt (PP_OUT+QQ_OUT); for t = 1: n if APPARENT_POWER_IN(t,1) < 1000, APPARENT_POWER_TEST(t,1) =1; else APPARENT_POWER_TEST(t,1) = 0; endif; %Ensure BOTH power in AND power out is under limit. if APPARENT_POWER_OUT(t,1) < 1000, APPARENT_POWER_TEST(t,2) =1; else APPARENT_POWER_TEST(t,2) = 0; endif; APPARENT_POWER_TEST(t,3) = APPARENT_POWER_TEST(t,2) + APPARENT_POWER_TEST(t,1); endfor %----------------------------------------------------------------------------------------------------------------% %PLOTTING THE RESULTS %----------------------------------------------------------------------------------------------------------------% %Set Viewpoint characteristics (i.e. Lat, Long, zoom) ZOOM LEVEL MUST BE BETWEEN 6 & 8 %V_P = [-33.75, 150.5, 8]; %Reasonable co-ordinates for zoom level 8 (Sydney level) %V_P = [-33, 149, 7]; %Reasonable co-ordinates for zoom level 7 (North Coast level) %V_P = [-33, 148, 6]; %Reasonable co-ordinates for zoom level 6 (NSW level) %----------------------------------------------------------------------------------------------------------------% %Get Image from Google Maps API LAT = num2str(V_P(1,1)); LONG = num2str(V_P(1,2)); ZOOM = num2str(V_P(1,3)); CENTER = strcat (LAT,",",LONG); mapbgf = urlread(strcat("http://maps.google.com/maps/api/staticmap?center=",CENTER,"&zoom=",ZOOM,"&size=640x640&sensor=false&style=feature:all",char(124),"visibility:off&style=feature:water",char(124),"saturation:-100",char(124),"lightness:-50",char(124),"visibility:on&style=feature:administrative.province",char(124),"element:geometry",char(124),"visibility:on&style=feature:landscape",char(124),"element:geometry",char(124),"lightness:100",char(124),"visibility:on")); bgfid=fopen("static.png","w"); fwrite(bgfid,mapbgf,"char"); fclose(bgfid); I = imread("static.png"); if P_N == 4, PLOT_GRID = 2; else PLOT_GRID = 1; endif subplot(PLOT_GRID,PLOT_GRID,PLOTS); imshow(I) hold on %----------------------------------------------------------------------------------------------------------------% %Zoom Level transform variables T = [ 6.05, 2.97, 1.452; %Transform matrix 6.9, 3.5, 1.75; 45.87, 90.9, 183.15; 53.145, 108.67, 221.23 ]; Z_n = V_P(1,3) - 5; %Zoom Number (1,2 or 3, where 1 = "zoom level 6", 2 = "zoom level 7", 3 = "zoom level 8") Y = -V_P (1,1) - T(1,Z_n); %translation in y X = V_P(1,2) - T(2,Z_n); %translation in x XD = T(3,Z_n); %dilation factor x YD = T(4,Z_n); %dilation factor y %----------------------------------------------------------------------------------------------------------------% %Color scheme %Ok lines- 4284d3 – 66,132,211, width 1 LINES = [11,97,164]/256; %Err lines - FFB540 – 255,181,64, width 2 ERRLINES = [255,191,0]/256; %Out lines – 000000 – 0,0,0 – width 1, could we make these dashed? OUTLINES = [255,49,0]/256; %Ok bus - 04346C – 4, 52,108 – Radius 1 - Dot BUS_COL = [3,62,107]/256; %Err bus - A66500 – 166,101,0 – Radius 2 – Concentric circles ERR_BUS_COL = [166,124,0]/256; %----------------------------------------------------------------------------------------------------------------% subplot(PLOT_GRID,PLOT_GRID,PLOTS); %Plot Branches for i = 1:n c = [sscanf(char(branch(5,i)), '%f,'),sscanf(char(branch(8,i)), '%f,')]; if ((APPARENT_POWER_IN(i,1)) == 0), plot (((c(2,:)-X)*XD),((Y+c(1,:))*-YD),"color",OUTLINES,"linewidth",2,"linestyle","--") elseif (abs (APPARENT_POWER_TEST(i,3)) < 2), plot (((c(2,:)-X)*XD),((Y+c(1,:))*-YD),"color",ERRLINES,"linewidth",2) else plot (((c(2,:)-X)*XD),((Y+c(1,:))*-YD),"color",LINES,"linewidth",1) endif endfor %Plot Buses for i = 1:B_n GL = sscanf(char(bus(12,i)), '%f,'); if results.bus(i,8) > 1.1 for k = 1:3 plot (((GL(2,1)-X)*XD),((Y+GL(1,1))*-YD), "o","color",ERR_BUS_COL, "markersize",(2*k)) endfor elseif results.bus(i,8) < 0.9 for k = 1:3 plot (((GL(2,1)-X)*XD),((Y+GL(1,1))*-YD), "o","color",ERR_BUS_COL,"markersize",(2*k)) endfor else for k = 1:3 plot (((GL(2,1)-X)*XD),((Y+GL(1,1))*-YD), "o","color",BUS_COL,"markersize",k) endfor endif endfor %----------------------------------------------------------------------------------------------------------------% %Plot Legend L_L = 3; % Line length in legend if P_N == 4, L_S = 20; else L_S = 10; endif; % Legend Spacing if P_N == 4, L_TS = 30; else L_TS = 15; endif; % Legend Title spacing %Legend images for k = 1:3 plot (L_P(1,1), L_P(1,2)+L_TS, "o","color",ERR_BUS_COL,"markersize",(2*k)) endfor for k = 1:3 plot (L_P(1,1), L_P(1,2)+L_TS+L_S, "o","color",BUS_COL,"markersize",k) endfor plot ([(L_P(1,1)-L_L),(L_P(1,1)+L_L)],[L_P(1,2)+L_TS+2*L_S,L_P(1,2)+L_TS+2*L_S],"color",OUTLINES,"linewidth",2) plot ([(L_P(1,1)-L_L),(L_P(1,1)+L_L)],[L_P(1,2)+L_TS+3*L_S,L_P(1,2)+L_TS+3*L_S],"color",ERRLINES,"linewidth",2) plot ([(L_P(1,1)-L_L),(L_P(1,1)+L_L)],[L_P(1,2)+L_TS+4*L_S,L_P(1,2)+L_TS+4*L_S],"color",LINES,"linewidth",1) %Legend text text (L_P(1,1), L_P(1,2), "LEGEND") text (L_P(1,1)+10, L_P(1,2)+L_TS, "Bus Voltage Outside Limits") text (L_P(1,1)+10, L_P(1,2)+L_TS+L_S, "Bus OK") text (L_P(1,1)+10, L_P(1,2)+L_TS+2*L_S, "Line Outage") text (L_P(1,1)+10, L_P(1,2)+L_TS+3*L_S, "Line Overload") text (L_P(1,1)+10, L_P(1,2)+L_TS+4*L_S, "Line Ok") %TITLE SCEN_TITLE = strcat ("Scenario",{" "},num2str(SN-1)); text (300,0, SCEN_TITLE) %----------------------------------------------------------------------------------------------------------------% endfor hold off 0.02 '0 MVA' '25 yr' '0 MW' 0 VoltageAngle ShuntConductance ReplacementPeriod New South Wales Australia Nation Nation State Region ReactivePowerOutput OperationCostAsPercentageOfFirstCost Octave MinimumVoltageMagnitude MinimumRealPowerOutput MaximumVoltageMagnitude MaximumRealPowerOutput MachineBaseMva LossZone Emissions 0 Has min cardinality DistributedEfficiency Calculating Demand for NGACs English NSW Greenhouse Gas Reduction Scheme GHG Emissions Existing Energy Calculating Demand for NGACs - GGAS Fact Sheet 2004 Year URL 1 1 Has max cardinality Title Subject Source NSW Government Organization Language Document Matthew Sullivan Moxy Knowledge Management Matthew Sullivan LastName FirstName Employer Expert Cite CitedBy InverseOf CentralizedEfficiency CapitalCost BusArea BaseVoltage Electricity Network